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1.0  Objective 


The  objective  of  this  work  is  to  identify  isotopic  ratios  suitable  for  analysis  via  mass 
spectrometry  that  distinguish  between  commercial  nuclear  reactor  fuel  cycles,  fuel  cycles  for 
weapons  grade  plutonium,  and  products  from  nuclear  weapons  explosions.  Methods  will  also  be 
detennined  to  distinguish  the  above  from  medical  and  industrial  radionuclide  sources.  Mass 
spectrometry  systems  will  be  identified  that  are  suitable  for  field  measurement  of  such  isotopes 
in  an  expedient  manner. 

2.0  Scope 

This  proposal  is  in  support  of  the  Basic  Research  Program  for  Combating  Weapons  of  Mass 
Destruction,  and  the  assessment  of  non-traditional  isotopic  ratios  by  mass  spectrometry  for 
analysis  of  nuclear  activities.  Isotopic  ratios  will  be  calculated  for  radionuclides  produced  in 
commercial  nuclear  reactor  fuel  cycles,  fuel  cycles  for  weapons  grade  plutonium,  and  nuclear 
weapons  explosions.  The  isotopic  ratios  that  best  identify  the  source  of  the  radionuclides  will  be 
selected.  Isotopes  with  a  combination  of  low  production  yields  and  low  mass  spectrometry 
detection  limits  will  be  removed  from  consideration.  Mass  spectrometry  techniques  will  be 
evaluated  as  a  function  of  their  ability  to  detect  and  qualify  the  radionuclides  of  concern.  A  mass 
spectrometry  system  design  will  be  identified  that  has  the  detection  sensitivity  necessary  for  the 
work  and  is  capable  of  field  operations. 

Option  Year  1(FY  2010)  specifically  focused  on  assessing  the  uncertainty  and  range  of  the 
results  obtained  in  the  first  two  years  of  this  project.  Similar  calculations  were  performed,  but 
utilizing  more  advanced  codes  (such  as  MCNPX)  and  different  databases.  Current  commercial 
and  future  proposed  fuel  cycles  were  modeled.  An  effort  will  also  be  conducted  to  develop  a 
test  data  set  to  extract  forensic  content  from  the  test  datasets.  Both  real  and  simulated  data  sets 
were  utilized.  Sample  preparation  methods  will  be  developed  for  measurements  in  Option  Year 
2  of  this  project. 

Option  Year  2  (FY  2011)  will  focus  on  the  quantification  capabilities  of  mass  spectrometry 
methods  for  the  elements  of  choice.  Forensic  identification  algorithms  and  a  software  tool  for 
forensic  analysis  will  also  be  developed  and  prototyped.  Mass  spectrometry  measurements  of 
the  radioactive  isotopes  of  forensic  interest  will  be  conducted. 


3.0  Background 

Isotopic  Ratios 

There  are  many  sources  for  radionuclides  in  our  environment.  These  include  naturalsources,  the 
commercial  nuclear  industry,  nuclear  weapons,  the  medical  industry,  and  other  sources.  Often 
times,  the  source  of  the  radionuclide  may  be  determined  through  just  identification  of  the 
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radionuclide.  If  radionuclides  are  produced  through  different  sources,  the  identification  of  the 
source  is  complex.  In  order  to  ascertain  a  specific  source  for  attribution,  radionuclide  ratios  are 
often  employed. 

Production  yields  of  radionuclides  from  fission  are  a  function  of  many  variables  including:  the 
fissile  material,  the  energy  spectrum  of  the  neutron  flux,  the  magnitude  of  the  neutron  flux,  and 
the  duration  of  the  irradiation.  As  a  result,  the  ratios  of  certain  radionuclides  are  highly 
dependent  of  these  variables  and  may  be  utilized  to  distinguish  between  radionuclides  produced 
from  nuclear  weapons,  medical  waste,  short  nuclear  fuel  cycles  (e.g.  239Pu  production  fuel 
cycles),  and  long  nuclear  fuel  cycles  (e.g.,  commercial  nuclear  fuel  cycles).  While  the  above  is 
easily  stated,  the  difficult  part  is  to  determine  which  radionuclide  ratios  should  be  utilized  for 
best  forensic  value. 

As  an  added  complication,  nuclear  debris  taken  for  forensic  analysis  often  does  not  come  directly 
from  the  source.  There  is  often  some  type  of  chemical  process  or  other  process  that  may  alter  the 
sample  composition.  Chemical  fractionation  issues  result  and  may  significantly  alter  ratios  of 
radionuclides  of  different  elements.  To  mitigate  this  problem,  it  is  best  to  examine  isotopic  ratios 
of  individual  elements  since  these  ratios  will  be  largely  unaltered  by  chemical  processes. 

Mass  Spectrometry 

Traditional  methods  for  radionuclide  detection  depend  upon  measuring  the  energy  released 
during  radioactive  decay.  Decay  counting  is  relatively  simple,  but  sample  prep  and  analysis  takes 
time  to  complete.  If  short-lived  radionuclides  have  already  decayed,  traditional  counting  can  be 
quite  slow.  Mass  spectrometry  (MS)  techniques  often  require  comparable  sample  prep  to  decay 
counting,  but  analysis  is  faster  since  MS  counts  atoms  rather  than  waiting  for  them  to  decay. 
Reducing  time  between  sample  collection  in  the  field  and  reliable  analytical  results  requires 
switching  to  MS. 

There  are  numerous  MS  techniques  capable  of  measuring  isotopic  ratios.  The  sample  size, 
detection  limit,  dynamic  range,  sample  prep  requirements,  and  ease  of  analysis  vary  widely 
among  the  techniques.  Some  techniques  have  very  simple  sample  prep,  requiring  only 
dissolution  in  acid  or  combustion  prior  to  analysis.  Others  require  extensive  preprocessing  that 
impedes  quick  turnaround.  In  practice,  the  selected  MS  technique  will  need  to  accurately 
measure  isotope  ratios  in  a  range  of  interest  as  quickly  as  possible.  To  insure  speedy  analysis,  the 
MS  technique  should  probably  be  sufficiently  robust  to  be  field  deployable  inside  a 
transportainer. 

Results  from  Curren  t  Work 

Extensive  work  was  conducted  during  the  first  two  years  of  this  project.  Nuclear  reactor  fuel 
cycles  were  modeled  utilizing  ORIGEN.  Fuels  cycles  from  pressurized  water  reactors  (PWR), 
boiling  water  reactors  (BWR),  and  Canadian  natural  deuterium  (CANDU)  reactors  were  all 
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evaluated.  Nuclear  weapons  were  modeled  by  utilizing  a  bare  sphere  (keff =1.0)  in  MCNPX 
utilizing  the  BURN  card.  The  production  of  every  fission  product,  activation  product,  and 
transuranic  was  recorded  and  entered  into  a  database. 

An  R  value  was  calculated  for  each  possible  isotopic  ratio.  This  is  a  metric  to  evaluate  the 
forensic  value.  R  values  greater  than  100  or  less  than  0.01  are  considered  good.  R  values  are 
calculated  as  shown  in  equation  1 . 


(1) 


Isotopic  ratios  were  then  prioritized  by  magnitude  of  the  ratio,  the  absence  of  possible 
interferences  in  field  monitoring,  and  the  mass  of  isotope  produced.  Table  1  shows  the  isotopic 
ratios  identified  as  the  best  for  distinguishing  between  an  unknown  reactor  type  and  a  known 
commercial  reactor  signature.  For  this  study  the  ratio  of  134Cs/135Cs  was  determined  to  be  the 
optimum  ratio. 


Table  1.  Top  forensic  indicators  to  differentiate  between  nuclear  weapons  and  commercial 
nuclear  reactors. 
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R  value  calculations  were  also  conducted  on  short  and  long  nuclear  fuel  cycles  in  PWR,  BWR, 
and  CANDU  reactors.  The  short  cycle  was  defined  as  one  that  produces  weapons  grade  Pu.  The 
long  fuel  cycle  was  defined  as  one  that  was  indicative  of  normal  commercial  nuclear  reactor 
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operation.  Similar  to  above,  the  ratios  were  then  prioritized  by  the  magnitude  of  the  R  value,  the 
absence  of  possible  interferences  in  field  monitoring,  and  the  mass  of  isotope  produced.  Table  2 
shows  the  results  for  the  best  forensic  indicators  to  determine  fuel  cycle  length.  The  ratio  of 
146Sm/151Sm  was  detennined  to  provide  the  best  forensic  value. 


Table2. .  Top  forensic  indicators  to  differentiate  between  commercial  nuclear  reactor  fuel 
cycle  length. 
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As  a  result  of  the  above  work,  Sm  and  Cs  were  identified  as  the  best  elements  to  focus  on  for 
nuclear  forensics  with  mass  spectrometry. 


4.0  Tasks/Scientific  Goals 
Option  Year  1  (FY  2011-2012) 

Task  4.1:  Utilize  MCNPX  to  calculate  isotopic  ratios  for  PWR,  BWR,  and  CANDU 
reactors.  Compare  results  with  those  from  ORIGEN  obtained  in  the  initial  Phase  of  this 
work. 
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The  Oak  Ridge  Isotope  Generation  and  Determination  code,  or  ORIGEN,  is  frequently  used  to 
detennine  nuclide  production  and  fuel  bum  up  in  a  nuclear  reactor.  This  work  sought  to  explore 
the  effectiveness  of  using  Monte-Carlo  simulations  to  perform  this  same  task.  Specifically,  the 
several  reactor  designs  were  run  using  MCNPX  to  calculate  fuel  burnup  and  fission  products,  as 
well  as  their  respective  radiogenic  daughters. 

This  work  compared  nuclide  production  using  ORIGEN  and  MCNPX  calculations  in 
three  reactor  types:  the  Westinghouse  17x17  pressurized  water  reactor,  the  General  Electric  8  x 
8  boiling  water  reactor,  and  the  CANDU-37  reactor.  Further,  a  simple  sensitivity  study  was 
performed  using  the  boiling  water  reactor  to  determine  what  effects  small  perturbations  had  on 
the  isotopic  ratios  calculated. 

For  this  task,  nuclide  production  during  fuel  irradiation  was  determined  using  pre-existing 
reactor  configurations  for  a  simple  boiling  water  reactor,  a  pressurized  water  reactor  and  a 
CANDU  reactor.  Nuclide  production  calculations  were  made  after  a  fuel  bumup  of  1  month  and 
18  months.  The  cycle  ratio,  Rc,  was  determined  for  a  number  of  nuclide  pairs  and  is  defined  by 
equation  1  above. 

Nuclide  production  was  then  detennined  using  probabilistic  Monte  Carlo  simulations,  as 
opposed  to  the  deterministic  ORIGEN  simulations  also  performed.  This  was  achieved  using  the 
BURN  card  in  MCNPX.  MCNP  input  decks  of  the  same  three  reactor  types  used  in  the 
ORIGEN  calculations:  the  CANDU-37,  the  GE  8  x8  BWR,  and  the  Westinghouse  PWR,  were 
used  in  this  work.  Cycle  ratios  were  then  calculated  using  the  output  of  the  MCNPX  burnup 
calculations  and  then  compared  to  the  cycle  ratios  detennined  using  the  ORIGEN  code. 

A  sensitivity  study  on  what  effects  small  changes  in  reactor  parameters  had  on  the  isotope  ratios 
was  performed.  The  sensitivity  measurements  were  conducted  using  MCNP  on  the  boiling  water 
reactor.  Changes  in  power  level,  bum  time,  and  initial  boron  concentrations,  were  investigated. 

In  order  to  minimize  discrepancies  between  the  MCNP  nuclide  production  output  and  the  results 
of  the  ORIGEN  runs,  efforts  were  made  to  ensure  that  the  reactor  parameters  used  for  both  the 
MCNPX  and  ORIGEN  input  files  were  the  same.  The  following  tables,  tables  3-8,  outline  the 
inputted  parameters  used  for  the  ORIGEN  portion  in  this  investigation.  The  fuel  mass,  type, 
enrichment,  and  also  the  moderator  density  were  extracted  from  the  respective  MCNPX  input 
decks.  The  parameters  stipulated  in  the  BURN  card  of  the  MCNPX  input  deck  also  matched  the 
bumup  and  reactor  power  levels  inputted  into  the  ORIGEN  run.  Figure  1  illustrates  the  MCNPX 
fuel  assembly  models  for  the  BWR,  PWR,  and  CANDU  models. 


Table  3:  Input  parameters  for  ORIGEN  simulation  of  BWR  following  1 -month  bumup 


FUEL  TYPE 

GE  8x8-4 

7 


U  INITIAL  MASS  (g) 

173,557 

ENRICHMENT 

3.23  w/o  235U 

BURNUP 

2952.34  MWd/MTU 

CYCLES 

1 

LIBRARIES 

1 

COOLING  TIME 

1  second 

MODERATOR  DENSITY 

0.6  g/cm3 

POWER 

17.08  MW 

AVERAGE  POWER 

98.411  MW/MTU 

Table  4:  Input  parameters  for  ORIGEN  simulation  of  BWR  following  18-month  bumup 


FUEL  TYPE 

GE  8x8-4 

U  INITIAL  MASS  (g) 

173,557 

ENRICHMENT 

3.23  w/o  235U 

BURNUP 

53142.20  MWd/MTU 

CYCLES 

1 

LIBRARIES 

1 

COOLING  TIME 

1  second 

MODERATOR  DENSITY 

0.6  g/cm3 

POWER 

17.08  MW 

AVERAGE  POWER 

98.411  MW/MTU 

Table  5:  Input  parameters  for  ORIGEN  simulation  of  PWR  following  1 -month  burnup 


FUEL  TYPE 

Westinghouse  17x17 

U  INITIAL  MASS  (g) 

450,030 

ENRICHMENT 

4.5  w/o  235U 

BURNUP 

1138.59  MWd/MTU 

CYCLES 

1 

LIBRARIES 

1 

COOLING  TIME 

1  second 

MODERATOR  DENSITY 

0.723  g/cm3 

POWER 

17.08  MW 

AVERAGE  POWER 

37.953  MW/MTU 

Table  6:  Input  parameters  for  ORIGEN  simulation  of  PWR  following  18-month  burnup 


FUEL  TYPE 

Westinghouse  17x17 

8 


U  INITIAL  MASS  (g) 

450,030 

ENRICHMENT 

4.5  w/o  235U 

BURNUP 

20494.63  MWd/MTU 

CYCLES 

1 

LIBRARIES 

1 

COOLING  TIME 

1  second 

MODERATOR  DENSITY 

0.723  g/cm3 

POWER 

17.08  MW 

AVERAGE  POWER 

37.953  MW/MTU 

Table  7:  Input  parameters  for  ORIGEN  simulation  of  CANDU  reactor  following  1 -month 
bumup 


FUEL  TYPE 

CANDU-37 

U  INITIAL  MASS  (g) 

19,832 

ENRICHMENT 

Natural 

BURNUP 

756.35  MWd/MTU 

CYCLES 

1 

LIBRARIES 

1 

COOLING  TIME 

1  second 

MODERATOR  DENSITY 

0.8121  g/cm3 

POWER 

0.5  MW 

AVERAGE  POWER 

2.88  MW/MTU 
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Table  8:  Input  parameters  for  ORIGEN  simulation  of  CANDU  reactor  following  1 -month 
bumup 


FUEL  TYPE 

CANDU-37 

U  INITIAL  MASS  (g) 

19,832 

ENRICHMENT 

Natural 

BURNUP 

1555.68  MWd/MTU 

CYCLES 

1 

LIBRARIES 

1 

COOLING  TIME 

1  second 

MODERATOR  DENSITY 

0.8121  g/cmJ 

POWER 

0.5  MW 

AVERAGE  POWER 

2.88  MW/MTU 
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PWR  -  Westinghouse  17x17 


Figure  1.  MCNPX  models  for  BWR,  PWR,  and  CANDU  fuel  assemblies. 


The  first  goal  of  this  work  was  to  determine  whether  R  values,  as  described  above,  have 
reasonable  agreement  for  the  radionuclides  of  interest  shown  in  Tables  1  and  2  of  this  report.  In 
addition,  a  sensitivity  study  was  conducted  to  investigate  how  the  R  values  change  over  a  range 
of  operational  conditions.  It  was  shown  that  the  R  values  compare  between  the  MCNPX  and  the 
ORIGEN  codes.  In  addition,  the  largest  operational  variable  that  appears  to  affect  the  R  values 
was  the  boron  concentration  in  the  moderator  of  the  PWR. 


Task  4.2:  Develop  fuel  cycle  models  that  predict  forensic  signatures  from  known  generic 
fuel  cycles 
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To  identify  isotopic  ratios  that  could  be  used  to  differentiate  short  (low  burnup)  from  long  (high 
bumup)  fuel  cycles,  ORIGEN-ARP  within  SCALE  6  was  used  to  model  several  fuel  types  and 
generate  expected  nuclide  compositions  that  would  arise  from  a  low  and  a  high  burnup  case. 
Typical  enrichments,  specific  power,  and  typical  bumup  values  for  commercial  applications 
were  gathered  from  available  sources.1'2  From  the  range  of  typical  burnup  values  upon 
discharge,  representative  values  for  a  low  and  a  high  burnup  case  were  estimated.  Additionally, 
some  simplifications  to  the  sample  parameters  were  made.  For  example,  for  the  specific  power, 
different  values  were  used  for  the  BWR  and  PWR  classes  rather  than  each  different  assembly 
type.  All  simulations  were  performed  using  “Express  Mode.”  The  sample  space  is  summarized 
in  Table  9. 

Table  9.  Summary  of  Samples  simulated  in  ORIGEN-ARP.  Note,  the  names  of  the  fuel  types 
follow  the  convention/format  of  ORIGEN-ARP 


Fuel  Type 

Reactor  Type 

Specific  Power 
[MW/MTU] 

Initial 

Enrichment 

[w/o] 

Burnup  Values 
[MWd/kgU] 

Abb8x8-1 

BWR 

23 

2.9 

1,27 

Atrium9-9 

BWR 

23 

2.9 

1,27 

Atrium  10-9 

BWR 

23 

2.9 

1,27 

CE14xl4 

PWR 

32 

2.8 

1,32 

CE16xl6 

PWR 

32 

2.8 

1,32 

GE7x7-0 

BWR 

23 

2.9 

1,27 

GE8x8-4 

BWR 

23 

2.9 

1,27 

GE9x9-7 

BWR 

23 

2.9 

1,27 

GE  10x1 0-8 

BWR 

23 

2.9 

1,27 

S14xl4 

PWR 

32 

2.8 

1,32 

Svea64-1 

BWR 

23 

2.85 

1,27 

Svea 100-0 

BWR 

23 

2.85 

1,27 

Vver440(3.6) 

PWR 

32 

3.60 

1,32 

Vver440(3.82) 

PWR 

32 

3.82 

1,32 

Vver440(4.25) 

PWR 

32 

4.25 

1,32 

Vver440(4.38) 

PWR 

32 

4.38 

1,32 

W14xl4 

PWR 

32 

2.8 

1,32 

W15xl5 

PWR 

32 

2.8 

1,32 

W17xl7 

PWR 

32 

2.8 

1,32 

1  Knief,  Ronald  Allen,  Nuclear  Engineering:  Theory  and  Technology  of  Commercial  Nuclear  Power,  2nd  Ed.,  American 
Nuclear  Society,  IL:  2008 

2  http://wp.ornl.gov/sci/scale/pubs/152495.pdf 
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After  simulation,  the  output  files,  which  each  contained  the  top  200  nuclides  (by  mass  at 
discharge)  were  each  processed  to  compute  all  isotopic  ratios  of  the  form 


7X  , 


I A ' 


X 


,  where  A  ^  A’.  Then,  these  ratios  as  defined  in  the  statement  of  work  were  compared 


between  the  low  and  high  bumup  cases  for  each  fuel  to  find  R-values  as  shown  in  equation  2: 

R, 


R  = 


—  n low  . 


f  D 

^ high 


(2) 


Lastly,  R-values  greater  than  100  or  less  that  0.01,  indicating  two  orders  of  magnitude  change  in 
the  isotopic  ratio  between  the  two  bumup  values,  were  sought.  Isotopic  ratios  with  such  R- 
values  are  identified  as  being  good  candidates  for  differentiating  a  low  and  high  bumup  fuel 
cycle.  These  identified  R-values  are  shown  in  Table  10. 


Task  4.3:  Inverse  calculations  to  assess  forensic  capabilities  from  test  data  sets 

This  task  was  comprised  of  three  major  steps:  simulate  fuel  cycles  and  assemble  nuclide 
inventories  into  a  database,  write  an  algorithm  to  compare  test  cases  against  entries  in  the 
database,  assigning  a  figure  of  merit  to  convey  their  similarity  and  determine  the  best-fit  entry  of 
the  database,  and  simulate  test  cases  to  test  the  identification  methodology’s  ability  to  identify  a 
known  fuel  cycle  that  best  matches  the  unknown. 


Creating  the  database 

To  mn  the  large  number  of  samples  required  for  the  assembly  of  the  desired  database,  a 
MATLAB  function  was  written  that  could  take  in  parameters  such  as  file  name,  fuel  type, 
enrichment,  cooling  time,  output  units,  etc.  and  write  an  ORIGEN-S  input  file  (.inp)  reflecting 
these  parameters.  Additionally,  a  script  was  written  to  call  this  function  for  every  point  in  the 
sample  space  and  write  the  required  batch  file  to  run  all  of  the  generated  input  files.  These 
functions  are  “makeOrigenlnp.m”  and  “makelnpsAndBat.m”  (see  Appendix).  These  files  were 
then  run  with  ORIGEN-S  within  SCALE  6.  All  samples  were  simulated  at  four  different  cooling 
times:  1  minute,  1  day,  30  days,  and  1  year.  Table  11  summarizes  the  sample  space  that  was 
simulated.  Note,  the  MATLAB  notation  for  an  array  of  values  is  used  for  denoting  the  burnup 
values:  600:690:13000  denotes  the  set  from  600  to  13000  in  steps  of  690.  The  step  sizes  were 
chosen  to  correspond  to  roughly  1  month  of  operation.  As  in  Task  4.2,  these  parameters  were 
estimated  from  available  literature  to  simulate  some  of  the  major  differences  between  these 
different  reactor/assembly  types,  while  keeping  the  simulations  simple. 
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Table  10.  Summary  of  identified  isotopic  ratios  useful  for  differentiation  between  low  and  high  burnup  cases.  Marked  are  ratios 
having  R-values  larger  than  100  or  less  than  0.01  for  a  particular  fuel  type  (continued  on  next  page). 
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Table  11.  Summary  of  samples  simulated  in  ORIGEN-ARP.  Note  the  names  of  the  fuel  types 
follow  the  convention/format  of  ORIGEN-ARP. 


Fuel  Type 

Reactor 

Type 

Specific  Power 
[MW/MTU] 

Initial 

Enrichment 

[w/o] 

Bumup  Values 
[MWd/kgU] 

Abb8x8-1 

BWR 

23 

2.9 

600:690:13000 

Atrium9-9 

BWR 

23 

2.9 

600:690:13000 

Atrium  10-9 

BWR 

23 

2.9 

600:690:13000 

CANDU28 

PHWR 

22 

0.711 

600:690: 130003 

CANDU37 

PHWR 

22 

0.711 

600:690:13000 

CE14xl4 

PWR 

32 

2.8 

600:950:18000 

CE16xl6 

PWR 

32 

2.8 

600:950:18000 

GE7x7-0 

BWR 

23 

2.9 

600:690:13000 

GE8x8-4 

BWR 

23 

2.9 

600:690:13000 

GE9x9-7 

BWR 

23 

2.9 

600:690:13000 

GE  10x1 0-8 

BWR 

23 

2.9 

600:690:13000 

S14xl4 

PWR 

32 

2.8 

600:950:18000 

Svea64-1 

BWR 

23 

2.85 

600:690:13000 

Svea 100-0 

BWR 

23 

2.85 

600:690:13000 

Vver440(3.6) 

PWR 

32 

3.60 

600:950:18000 

Vver440(3.82) 

PWR 

32 

3.82 

600:950:18000 

Vver440(4.25) 

PWR 

32 

4.25 

600:950:18000 

Vver440(4.38) 

PWR 

32 

4.38 

600:950:18000 

VverlOOO 

PWR 

32 

2.8 

600:950:18000 

W14xl4 

PWR 

32 

2.8 

600:950:18000 

W15xl5 

PWR 

32 

2.8 

600:950:18000 

W17xl7 

PWR 

32 

2.8 

600:950:18000 

To  fonn  the  database,  a  master  list  of  nuclides  present  in  the  1628  samples  was  fonned.  This  list 
includes  308  nuclides.  A  308  x  1628  matrix,  D,  was  then  created,  where  each  row  corresponds 
to  a  nuclide  in  the  master  list  and  each  column  a  sample.  The  matrix  was  populated  as  follows: 

D(/  =  Mass  of  Nuclide  i  in  Sample  j  (3) 

where  the  nuclide  index  i  is  taken  from  the  master  list.  If  sample  j  did  not  have  nuclide  i  in  its 
top-200  list,  a  zero  was  placed  in  that  entry  of  D  (see  “makeDTRAdatabase.m”  in  Appendix). 

3  These  burnup  values  are  too  high  for  CANDU  reactors,  but  for  simplicity,  the  same  values  as  used  for  other 
reactors  with  a  similar  power  density  were  used.  The  effects  of  heavy-water  moderation,  cross  sections,  etc.  are 
still  simulated. 
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Writing  an  identification  algorithm 


To  identify  an  unknown  fuel  cycle  from  its  top-200  nuclide  list,  an  algorithm  was  written  that 
sorts  the  unknown  cycle’s  list  to  agree  with  the  master  nuclide  list  described  above  (inserting 
zeros  as  appropriate),  then  computes  a  figure  of  merit  (FOM)  to  describe  the  similarity  between 
the  resulting  column  vector  and  the  1628  columns  of  D.  This  was  done  by  three  methods: 


FOM1^YJl0ffxn-xny 
FOM2  =  Fn°=\ \xn-xn\ 
FOM3  =  £308  fnzfA 


(4) 

(5) 

(6) 


Here  xn  is  the  nth  entry  of  the  database  matrix  for  a  particular  sample  and  xn  is  the  nth  entry  of 
the  unknown  sorted  nuclide  list.  These  three  FOM’s  each  have  a  motivation:  1  is  standard  least- 
squares  approach  that  is  commonly  used  in  many  “best  fit”  applications,  2  adjusts  the  approach 
of  1,  which  may  over  penalize  a  dissimilar  entry,  and  three  adjusts  this  second  approach  to 
examine  the  relative  difference  in  each  entry  in  order  to  more  fairly  weight  the  contributions  to 
the  FOM  from  low  concentration  nuclides.  In  practice,  an  FOM  is  generated  for  each  of  the 
1628  columns  of  D,  and  the  fuel  cycle  (fuel  type,  cooling  time,  and  burnup)  with  the  smallest 
entry  in  the  resulting  vector  is  found  to  be  the  best  fit.  Additionally,  an  estimate  of  the  burnup  of 
the  unknown  is  made  by  doing  a  linear  interpolation  between  two  samples  nearest  to  the  best-fit 
entry  in  D  (which  was  organized  to  place  column  vectors  of  nuclide  inventories  of  simulations  of 
the  same  reactor  and  cooling  time  next  to  each  other  in  order  of  bumup).  The  FOM  was  used  as 
the  notion  of  distance  (see  “IDfuelCycle.m”  in  Appendix). 

To  verify  the  quality  of  the  code,  samples  from  the  database  were  put  into  the  identification 
algorithm  in  order  to  verify  that  the  correct  column  was  returned  as  the  best  fit  and  the 
corresponding  FOM  value  was  0  as  would  be  expected.  This  was  done  using 
“checkDatabase.m.”  In  this  process,  it  was  realized  that  the  third  FOM  was  ill-defined  (the 
inclusion  of  zeros  in  the  database  causes  the  FOM  column  vectors  to  have  undefined  entries). 
Rather  than  change  the  database  formulation,  this  third  FOM  option  was  not  further  used. 


Testing  the  method 
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To  test  the  method,  a  test  set  of  ORIGEN  simulations  was  performed  using  the  aforementioned 
MATLAB  routines.  These  test  cases  featured  burnup  values  that  were  not  exactly  represented  in 
the  database,  i.e.,  falling  in  between  bumup  values  of  sample  populating  the  database.  Similar 
perturbations  were  also  made  on  the  enrichment  of  the  fuel,  cooling  time  of  the  fuel,  and  in  the 
various  combinations  of  these  three  variables.  Table  12  summarizes.  These  samples  were  fed 
into  the  identification  function  for  comparison  against  the  database  matrix  and  the  results 
analyzed. 

Table  12.  Summary  of  the  test  cases  used  for  testing  the  identification  algorithm. 


Fuel 

Enrichment  [w/o] 

Cooling  Time 

Burnup  [MWd/MTU] 

CANDU28 

0.71 

lm,  7d,  30d,  la 

1400,  5000,  11000 

CANDU28 

0.71 

3m,  9d,  2a 

6120,  5000 

CE  16x16 

2.8 

lm,  7d,  30d,  la 

1700,  8700,  17000 

CE  16x16 

2.8 

3m,  9d,  2a 

8700,9150 

CE  16x16 

3.1 

7d,  9d 

8700,9150 

GE  7x7-0 

2.9 

lm,  7d,  30d,  la 

2000,  7200,  10800 

GE  7x7-0 

2.9 

3m,  9d,  2a 

7200,  8880 

GE  7x7-0 

3.2 

7d,  9d 

7200,  8880 

The  function  “runTestCases.m”  performed  these  tests  and  recorded  the  results.  The  names  of  the 
samples  were  parsed  to  find  the  parameters  used  in  the  simulation,  and  these  parameters  (bumup, 
enrichment,  cooling  time)  were  examined  to  determine  which  were  perturbed  versus  calibrated 
values.  The  function  then  passed  one  of  the  test  cases  into  the  identification  algorithm.  Relative 
errors  in  the  estimation  of  burnup,  both  from  the  best-fit  entry  of  the  database  matrix  and  from 
the  linear  interpolation  method,  were  calculated  for  each  test  case.  Lastly,  successes  in  the 
estimation  of  the  fuel  assembly  type  and  cooling  time  were  recorded.  Successes  in  the 
estimation  of  the  assembly  type  were  defined  as  getting  the  type  exactly  correct,  i.e.,  guessing 
that  a  test  sample  was  a  CE  16x16  rather  than  CE  14x14  was  not  considered  a  success.  For  a 
cooling  time  estimation  to  be  declared  successful,  it  had  to  be  the  closest  possible  guess  (of 
which  the  algorithm  can  only  make  4:  1  minute,  1  day,  30  days,  and  1  year)  to  the  actual  cooling 
time  of  the  unknown. 

Figure  1  shows  the  values  of  FOMi  and  FOIVF  for  different  entries  in  the  database,  generated  for 
one  of  the  test  cases.  In  this  example,  the  simulated  fuel  assembly  was  a  CE  16x16  type,  with  an 
initial  enrichment  of  2.8  w/o,  a  burnup  value  of  1700  MWd/MTU,  and  a  cooling  time  of  one 
year.  The  algorithm  is  able  to  correctly  identify  the  type  of  reactor  (see  Figure  2. a)  and  then  the 
nearest  burnup  value  represented  in  the  database  (see  Figure  2.b).  Both  of  these  identifications 
are  indicated  by  a  minimum  in  the  FOM  values.  Note,  while  these  FOMs  are  discrete  functions, 
they  are  displayed  with  a  line  for  ease  of  reading  the  figures.  The  interpolated  burnup  value 
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produced  by  the  second  FOM  produces  a  better  estimate  of  the  burnup  of  the  test  case  versus  the 
best-fit  entry  from  the  database,  with  relative  errors  of  3.24%  and  8.82%,  respectively. 

These  and  the  rest  of  the  results  for  test  cases  that  only  altered  the  bumup  are  given  in  Table  12. 
The  complete  results,  as  described  above,  became  very  difficult  to  examine,  and  further  work  is 
needed  to  fully  understand  underlying  trends,  sensitivities  of  the  method  to  the  different  sources 
of  variance  introduced  in  this  study,  and  to  alter  the  method  to  better  deal  with  these  additional 
sources  of  variance.  The  complete  results  may  be  found  in  Tables  A.l  and  A.2  within  the 
Appendix. 
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Figure  2.  Values  of  FOM  for  entries  in  the  database,  l.a:  (Top)  The  area  around  the  minimum 
(black  line)  identifies  the  reactor  type  and  cooling  time,  l.b:  (Bottom)  Shows  the  FOM  values 
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for  entries  in  the  database  corresponding  to  CE  16x16  with  1  year  cooling  time  and  minimum 
vale  for  the  two  FOM's,  identifying  the  best  guess  for  the  bumup  of  the  test  case. 
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Table  7.  Shows  the  results  of  test  cases  where  only  burnup  was  altered.  The  results  shown  here  were  obtained  with  the  least  squares 
value  figure  of  merit.  The  case  corresponding  to  Figure  1  is  highlighted  and  shows  the  same  test  case  results  as  obtained  using  the 
absolute  value  figure  of  merit. 


Reactor 

Type  Guess 

Simulated 

Reactor 

Success 

Best-Fit 

Burnup 

[MWd/MTU] 

Error 

[%] 

Interpolated 

Burnup 

[MWd/MTU] 

Error 

[%] 

Simulated 

Burnup 

[MWd/MTU] 

Cooling 

Time 

Guess 

Simulated 

Cooling 

Time 

Success 

CANDU28 

CANDU 

28 

X 

1290 

7.86 

1552 

10.86 

1400 

la 

la 

X 

CANDU  37 

CANDU 

28 

1290 

7.86 

1522 

8.71 

1400 

7d 

lm 

CANDU  28 

CANDU 

28 

X 

1290 

7.86 

1547 

10.50 

1400 

30d 

7d 

CANDU  28 

CANDU 

28 

X 

1290 

7.86 

1552 

10.86 

1400 

30d 

30d 

X 

CANDU  37 

CANDU 

28 

4740 

5.20 

5318 

6.36 

5000 

la 

la 

X 

CANDU  37 

CANDU 

28 

4740 

5.20 

5301 

6.02 

5000 

7d 

lm 

CANDU  37 

CANDU 

28 

4740 

5.20 

5313 

6.26 

5000 

30d 

7d 

CANDU37 

CANDU 

28 

4740 

5.20 

5318 

6.36 

5000 

30d 

30d 

X 

CANDU  28 

CANDU 

28 

X 

10950 

0.45 

11046 

0.42 

11000 

la 

la 

X 

CANDU  28 

CANDU 

28 

X 

10950 

0.45 

11047 

0.43 

11000 

lm 

lm 

X 

CANDU28 

CANDU 

28 

X 

10950 

0.45 

11047 

0.43 

11000 

7d 

7d 

X 

CANDU  28 

CANDU 

28 

X 

10950 

0.45 

11047 

0.43 

11000 

30d 

30d 

X 
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CE  16x16 

CE  16x16 

X 

1550 

8.82 

1909 

12.29 

1700 

la 

la 

X 

CE  16x16 

CE  16x16 

X 

1550 

8.82 

1889 

11.12 

1700 

7d 

lm 

CE  16x16 

CE  16x16 

X 

1550 

8.82 

1905 

12.06 

1700 

30d 

7d 

CE  16x16 

CE  16x16 

X 

1550 

8.82 

1909 

12.29 

1700 

30d 

30d 

X 

VVER  1000 

CE  16x16 

9150 

5.17 

8307 

4.52 

8700 

la 

la 

X 

VVER  1000 

CE  16x16 

9150 

5.17 

8317 

4.40 

8700 

lm 

lm 

X 

VVER  1000 

CE  16x16 

9150 

5.17 

8334 

4.21 

8700 

lm 

7d 

VVER  1000 

CE  16x16 

9150 

5.17 

8309 

4.49 

8700 

30d 

30d 

X 

W  14x14 

CE  16x16 

16750 

1.47 

17198 

1.16 

17000 

la 

la 

X 

W  14x14 

CE  16x16 

16750 

1.47 

17206 

1.21 

17000 

lm 

lm 

X 

W  14x14 

CE  16x16 

16750 

1.47 

17196 

1.15 

17000 

30d 

7d 

W  14x14 

CE  16x16 

16750 

1.47 

17201 

1.18 

17000 

30d 

30d 

X 

GE  7x7-0 

GE  7x7-0 

X 

1980 

1.00 

2062 

3.10 

2000 

la 

la 

X 

GE  7x7-0 

GE  7x7-0 

X 

1980 

1.00 

2054 

2.70 

2000 

lm 

lm 

X 

GE  7x7-0 

GE  7x7-0 

X 

1980 

1.00 

2060 

3.00 

2000 

7d 

7d 

X 

GE  7x7-0 

GE  7x7-0 

X 

1980 

1.00 

2061 

3.05 

2000 

30d 

30d 

X 

GE  8x8-4 

GE  7x7-0 

7500 

4.17 

6922 

3.86 

7200 

la 

la 

X 

GE  8x8-4 

GE  7x7-0 

7500 

4.17 

6913 

3.99 

7200 

lm 

lm 

X 

GE  8x8-4 

GE  7x7-0 

7500 

4.17 

6932 

3.72 

7200 

lm 

7d 

GE  8x8-4 

GE  7x7-0 

7500 

4.17 

6938 

3.64 

7200 

lm 

30d 

GE  7x7-0 

GE  7x7-0 

X 

10950 

1.39 

10666 

1.24 

10800 

la 

la 

X 

GE  7x7-0 

GE  7x7-0 

X 

10950 

1.39 

10658 

1.31 

10800 

lm 

lm 

X 

GE  7x7-0 

GE  7x7-0 

X 

10950 

1.39 

10662 

1.28 

10800 

7d 

7d 

X 

GE  7x7-0 

GE  7x7-0 

X 

10950 

1.39 

10669 

1.21 

10800 

7d 

30d 
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The  second  figure  of  merit  (absolute  value)  was  more  successful  in  identifying  the  reactor  type, 
with  27  successes  versus  19  successful  identifications  made  by  the  first  FOM,  while  the  first 
FOM  performed  slightly  better  at  predicting  cooling  time,  the  first  and  second  FOM’s  yielding 
25  and  24  successful  cooling  time  guesses  respectively.  The  average  relative  errors  in  predicting 
the  burnup  of  the  test  cases  using  the  best-fit  entry  from  the  database  was  similar  for  the  two 
FOM’s,  the  least  squares  averaging  3.95%  and  the  absolute  value  averaging  3.99%  relative  error. 
The  error  in  this  burnup  prediction  is  a  function  of  the  difference  between  the  bumup  values 
simulated  for  the  formation  of  the  database.  Unsurprisingly,  the  absolute  value  FOM  performed 
much  better  predicting  the  bumup  of  the  test  cases  via  linear  interpolation,  yielding  an  average 
relative  error  in  prediction  of  1.38%  versus  the  4.72%  average  produced  using  the  first  figure  of 
merit.  Thus,  overall  the  second  figure  of  merit  performed  better. 

As  the  full  set  of  test  cases  and  results  becomes  much  more  complex  to  analyze,  further  work  is 
suggested  to  more  fully  understand  the  performance  of  the  method  outlined  in  this  work  and 
develop  ways  to  overcome  its  shortcomings.  Possible  solutions  could  include  using  different 
figures  of  merit  and/or  combinations  of  different  FOM’s  for  the  various  parts  of  the  analysis. 
Optimization  of  the  database,  e.g.,  including  more  simulations  or  concentrating  database  entries 
at  points  in  the  sample  space  (reactor  type,  burnup  values,  enrichment,  etc.)  where  the  output  is 
most  sensitive  to  changes  in  these  input  parameters,  could  also  improve  the  ability  of  the  method 
to  make  accurate  predictions.  Similarly,  the  sensitivity  of  each  nuclide  in  the  output  files  relative 
to  these  input  parameters  could  be  studied  and  this  knowledge  used  to  optimize  the  method. 
Lastly,  such  a  sensitivity  study  could  also  be  used  to  make  estimates  of  the  uncertainty  in  the 
various  predictions  made. 

In  summation,  the  present  work  demonstrates  the  efficacy  of  a  simple  method  of  analyzing  more 
of  the  data  available  from  measurements  of  spent  nuclear  fuel  for  the  purposes  of  fuel  cycle 
identification.  Further  work  is  warranted  to  develop  such  methods  as  a  complement  to  currently- 
employed  methods  that  rely  on  data  reduction  via  the  use  of  isotopic  ratios. 


Task  4.4:  Develop  chemistry  for  analysis  of  samarium. 

The  chemical  form  of  samarium  required  for  analysis  varies  for  different  mass  spectrometry 
techniques.  Inductively  Coupled  Plasma  Mass  Spectrometry  (ICPMS)  requires  a  liquid  sample 
for  introduction  into  the  plasma.  Accelerator  mass  spectrometry  requires  solid  targets  for 
negative  ion  production.  We  previously  examined  ion  production  from  samarium  oxide  (Sn^Cfi) 
and  samarium  fluoride  (SmFs)  in  the  LLNL  ion  source.  The  fluoride  produced  an  order  of 
magnitude  more  negative  ions  and  was  selected  as  the  best  option  for  samarium  AMS  sample 
material.  Samarium  is  readily  soluble  in  nitric  acid  or  hydrochloric  acid  and  commercial 
elemental  liquid  standards  contain  nitric  acid.  Liquid  Sm  standards  are  generally  produced  by 
dissolving  high  purity  solid  Sm203  in  high  purity  nitric  acid. 
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Production  of  SmF3  for  AMS  samples  was  investigated  starting  with  a  Sm  liquid  standard.  This 
material  does  not  have  elemental  interferences  and  is  a  good  starting  material  for  developing 
sample  chemistry.  The  LLNL  AMS  facility  uses  the  fluorides  of  Ca  and  Sr  currently,  so  we  have 
lab  facilities  suitable  for  working  with  hydrofluoric  acid.  The  general  procedure  for  making 
SmF3  is  as  follows:  Start  with  Sm  solution,  e.g.  4  mL  solution  containing  2  mg  Sm  in  5-10% 
nitric  acid.  Add  concentrated  ammonium  hydroxide  but  keep  solution  acidic.  Add  3  mL  40% 
hydrofluoric  acid  and  allow  samarium  fluoride  to  precipitate  overnight.  Centrifuge  to  form  a 
pellet  and  remove  the  supernatant.  Add  1  mL  of  dilute  hydrofluoric  acid,  then  transfer  samarium 
fluoride  and  rinse  to  a  new  2  mL  centrifuge  tube.  Centrifuge  again  to  form  a  pellet,  remove  the 
supernatant,  and  dry  overnight  on  a  heating  block.  Place  dry  SmF3  in  a  watertight  vial  until 
ready  to  use. 

Samarium  fluoride  can  be  precipitated  by  adding  a  soluble  fluoride  salt  such  as  ammonium 
fluoride  rather  than  hydrofluoric  acid.  Ammonium  fluoride  has  traces  of  hydrofluoric  acid  , 
however,  so  the  same  safety  procedures  are  used.  Hydrofluoric  acid  has  been  shown  to  work 
better  when  solutions  are  less  pure  (e.g.,  Ca),  however,  so  we  decided  to  stick  with  it. 


Task  4.5:  Conduct  presentations/meetings  at  times  and  places  specified  in  the  contract 
schedule. 

A  presentation  was  delivered  on  July  22,  2011  in  Springfield,  VA  at  the  DTRA  program  review. 


Task:  4.6:  Write  Option  Year  1  report 

It  was  written  and  submitted  by  the  September  1,  201 1  deadline. 


Option  Year  2  (FY  2011) 

Task  4.7  Validate  mass  spectrometry  detection  limits  for  identified  elements. 

Measurements  were  conducted  at  LLNL  utilizing  both  AMS  and  ICPMS  to  evaluate  the 
detection  capabilities  for  samarium  isotopes.  The  Sm203  sample  used  for  irradiation  appears  to 
be  extraordinarily  pure  samarium  for  a  commercial  material  as  analyzed  by  ICPMS.  No 
europium  or  15 ’Sm  was  identified  in  it.  Scans  on  a  pulse  counting  detector  over  mass  151  with 
simultaneous  collection  of  the  150Sm,  152Sm  and  154Sm  ion  beams  in  Faraday  cups  showed  only  a 
relatively  smooth  scattered-ion  tail  coming  from  152Sm.  No  indication  of  a  peak  could  be 
discerned  above  this  baseline,  which  was  about  50  cps  at  mass  151  for  a  4.534xl07  cps  ion  beam 
from  152Sm.  This  gives  an  abundance  sensitivity  at  mass  151  with  respect  to  152  of 
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approximately  1.1  ppm.  The  150Sm  ion  beam  was  running  at  about  1.257  xl07cps,  and  the  semi- 
quantitative  value  for  150Sm/152Sm  equals  0.277,  which  is  approximately  that  of  nonnal  Sm.  The 
best  estimate  of  a  count  rate  that  would  be  discernable  above  this  baseline,  the  critical  level,  is  25 
cps,  and  this  gives  a  result  for  151Sm/150Sm  <  2x1 0~6.  Hence,  the  detection  limit  for  ICPMS  of 
151Sm/150Sm  is  2x1  O'6. 

AMS  measurements  on  Sm  show  that  limitations  are  dominated  by  an  isobar  of  151Eu  and  a 
scattered  tail  of  stable  150Sm  isotope  that  cannot  be  distinguished  with  the  high  count  rates 
encountered.  The  l51Eu  interference  may  be  addressed  through  improved  chemistry.  The 
concentration  of  Eu  would  need  to  be  reduced  from  ~lxl0"6  to  lxlCT10  or  preferably,  lxlO'12. 

The  Sm  stable  isotope  interference  will  be  a  function  of  stable  Sm  in  the  sample.  For  pure 
nuclear  fuel  samples,  the  stable  Sm  isotope  interference  may  be  manageable.  For  environmental 
samples,  the  magnitude  of  the  interference  largely  depends  on  the  sample  makeup.  Sm  is  a 
lanthanide  with  an  average  crustal  composition  of  8  ppm.  As  a  result,  natural  Sm  interference 
with  anthropogenic  Sm  isotopes  could  cause  issues  depending  on  sample  acquisition. 

More  discussion  may  be  found  under  the  discussion  of  Task  4.9. 

Task:  4.8  Develop  algorithms  for  forensic  determinations  utilizing  mass  spectrometry. 

While  isotopic  ratios  have  been  used  effectively  until  present,  a  great  deal  of  information  available  to  the 
forensic  analyst  goes  unused  by  these  univariate  approaches.  Multivariate  analysis  has  the  ability  to 
overcome  difficulties  that  traditionally  plague  univariate  approaches,  such  as  poor  calibration  data  [6]. 
For  instance,  multivariate  data  analysis  allows  for  the  determination  of  multiple  analytes  with  a  single 
multivariate  measurement  and  the  detection  of  unexpected  components  [7].  Multi-way  analysis  can 
additionally  allow  for  adequate  analysis  of  analytes  in  the  presence  of  unexpected  components  and 
signals.  These  approaches  have  been  widely  applied  in  chemometrics  and  there  is  a  wealth  of  literature  on 
the  subject  [8,  9],  Additionally,  the  new  emphasis  on  speed  of  analysis  requires  as  much  useful 
information  to  be  gathered  from  every  measurement.  For  these  reasons,  we  believe  that  investigation  into 
the  application  of  multivariate  data  analysis  techniques  to  nuclear  forensic  applications  is  warranted. 

This  work  is  an  initial  study  to  validate  our  claim  that  multivariate  analysis  could  be  successfully  applied 
to  nuclear  forensic  applications.  We  do  this  by  applying  nearest  neighbor  and  ridge  regression  algorithms 
to  predict  the  burnup  and  reactor  type  of  spent  nuclear  fuel  using  simulated  data. 

Predictions  of  reactor  type  and  burnup  were  made  by  comparing  the  radionuclide  composition,  herein 
referred  to  as  the  nuclide  vector,  of  an  unknown  test  case  to  a  collection  of  vectors  corresponding  to 
known  reactor  type  and  burnup  value,  referred  to  as  the  training  data.  All  the  data  was  simulated,  and  the 
preparation  of  this  data  is  described  in  Section  3.1. 

The  training  data  was  loaded  into  a  matrix  X  G  MiVxp,  where  N  is  the  number  of  observations  and  p  is  the 
number  of  nuclides  whose  compositions  are  tracked. 
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Three  methods  were  used  to  compare  the  test  case  to  the  training  data: 

1 .  Nearest  Neighbor  with  L!  norm 

2.  Nearest  Neighbor  with  L2  norm 

3.  Ridge  Regression 

The  first  two  methods  utilize  the  nearest  neighbor  approach,  where  a  measure  of  dissimilarity/distance  is 
calculated  for  each  of  the  nuclide  vectors  in  the  training  data.  The  minimum  in  the  distance  measure  is 
found,  and  the  responsible  member  of  the  training  data  is  declared  the  best  fit.  The  reactor  type  and 
burnup  of  this  training  observation  is  taken  to  be  the  predictions  for  the  test  case. 

The  distance  measures  for  methods  one  and  two  are  given  as 

=  \\x-^\\\  =  j2\{xj-nf)\ 

j= i 

D&  =  \\x-^\\l  =  '£(xj-$))2 

.7  —  1 

where  x  is  the  test  case  nuclide  vector  and  //"  is  the  ith  column  of  the  training  data,  X. 

To  improve  the  fidelity  of  the  burnup  prediction  using  these  two  methods,  after  the  best- fit  member  of  the 
training  data  was  found,  a  linear  interpolation  was  performed  between  training  points  adjacent  to  the  best- 
fit  training  sample. 

Ridge  Regression  was  also  used  to  make  burnup  predictions.  This  is  a  standard  linear  regression  with  a 
penalty  placed  on  the  regression  vector,  /?  6  jr(JVx(p+1)  ),  to  discourage  overfitting3.  Given  the  true 
burnup  of  each  training  case,  y,  minimizing  the  ridge  regression  objective  gives  the  regression  vector. 

mmjlly—  [1  X]pf2  +  (9) 

where  1  is  column  vector  of  ones. 

Ridge  regression  was  chosen  because  it  has  a  closed  form  solution  and  penalizing  the  regression  vector 
makes  it  well  suited  for  applications  with  correlated  variables,  which  naturally  arises  in  the  analysis  of 
radionuclides  due  to  radioactive  decay  [9],  The  ridge  regression  model  was  trained  using  the  data  without 
added  error  (see  Section  3.2)  using  a  validation  set  that  was  removed  from  the  training  data. 

The  average  percent  error  in  burnup  prediction  was  recorded  as  a  measure  of  the  predictive  power  for 
each  combination  of  nuclide  subset  and  prediction  technique. 

Training  and  testing  data  was  generated  using  ORIGEN-ARP  5.1.01  included  in  SCALE  6.  The  top  200 
nuclides  by  mass  at  discharge  were  calculated  for  eleven  PWR  fuel  types,  nine  BWR  fuel  types,  and  two 
PWHR  fuel  types  using  the  cross  sections  included  in  ORIGEN-ARP.  Decay  calculations  were  performed 
for  four  different  cooling  times  after  discharge:  one  minute,  seven  days,  thirty  days,  and  one  year. 


(7) 

(8) 
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For  each  combination  of  cooling  time  and  fuel  type,  calculations  were  performed  for  a  wide  range  of 
burnup  values,  corresponding  to  approximately  one  month  to  eighteen  months  of  operation  in  one-month 
increments.  Power  densities  of  32,  23,  and  22  MW/MTU  were  used  for  the  simulations  of  the  PWRs, 
BWRs,  and  PHWRs,  respectively.  Burnup  values  used  for  PWR  simulations  were  600  to  17,700 
MWd/MTU  in  steps  of  950  MWd/MTU,  and  values  of  600  to  12,330  in  steps  of  690  MW d/MTU  were 
used  for  simulations  of  BWRs  and  PHWRs.  This  resulted  in  eighteen  burnup  values  for  the  BWR  and 
PHWR  reactors  and  nineteen  values  for  the  PWR  data.  The  basis,  enrichment,  and  availability  were  kept 
constant  between  calculations,  except  in  the  case  of  the  VVER  and  CANDU  reactors,  where  enrichment 
is  a  key  design  difference.  These  values  and  average  initial  enrichment  values  were  approximated  from 
the  literature  [10-12],  Table  8  summarizes  the  training  data  simulations. 

Table  8.  Parameters  used  in  ORIGEN-ARP  to  generate  training  data.  Fuel  types  are  referred  to  by  their 
identifiers  used  in  ORIGEN-ARP. 


Fuel 

Reactor 

Enrichment  [w/o] 

CE14xl4 

PWR 

2.8 

CE16xl6 

PWR 

2.8 

W14xl4 

PWR 

2.8 

W15xl5 

PWR 

2.8 

W17xl7 

PWR 

2.8 

S14xl4 

PWR 

2.8 

VVER440(3.60) 

PWR 

3.60 

VVER440(3.82) 

PWR 

3.82 

VVER440(4.25) 

PWR 

4.25 

VVER440(4.38) 

PWR 

4.38 

VVER1000 

PWR 

2.8 

GE7 

BWR 

2.9 

GE8 

BWR 

2.9 

GE9 

BWR 

2.9 

GE10 

BWR 

2.9 

Abb8 

BWR 

2.9 

Atrium9 

BWR 

2.9 

SVEA64 

BWR 

2.9 

SVEA100 

BWR 

2.9 

CANDU28 

PHWR 

0.711 

CAND37 

PHWR 

0.711 

Test  cases  were  generated  in  the  same  way.  Three  reactors  were  used,  and  initial  enrichments  and  cooling 
times  that  were  and  were  not  represented  in  the  training  data  were  simulated.  Table  9  summarizes  the  test 
cases.  In  total,  there  were  1,628  samples  in  the  training  set  and  62  in  the  testing  set. 
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Table  9.  Parameters  for  the  generation  of  test  case  data  in  ORIGEN-ARP.  Fuel  types  are  referred  to  by 
their  identifiers  used  in  ORIGEN-ARP.  Power  densities  are  the  same  as  the  training  data. 


Fuel 

Reactor 

Enrichment  [w/o] 

Cooling  Time 

Burnup  [MWd/MTU] 

CANDU28 

PHWR 

0.711 

{lm,  7d,  30d,  la} 

{1400,5000,11000} 

CANDU28 

PHWR 

0.711 

{3m,  9d,  2a} 

{5000,  6120} 

CE16xl6 

PWR 

2.8 

{lm,7d,  30d,la} 

{1700,  8700,  17000} 

CE16xl6 

PWR 

2.8 

{3m,  9d,  2a} 

{8700,  9150} 

CE16xl6 

PWR 

3.1 

{7d,  9d} 

{8700,  9150} 

GE7 

BWR 

2.9 

{lm,  7d,  30d,  la} 

{2000,  7200, 10800} 

GE7 

BWR 

2.9 

{3m,  9d,  2a} 

{7200,  8800} 

GE7 

BWR 

3.2 

{7d,  9d} 

{7200,  8800} 

The  top  200  nuclides  by  mass  at  discharge  were  recorded  for  each  case.  A  master  list  of  all  the  nuclides 
present  in  the  nuclide  vectors  within  the  training  data  were  pooled,  and  each  of  the  nuclide  vectors  in  the 
training  and  test  set  were  sorted  to  follow  this  order.  This  yielded  a  total  of  308  unique  radionuclides.  If  a 
nuclide  in  the  master  list  was  not  in  a  sample’s  top  200,  a  zero  was  used  in  that  position  of  the  formatted 
nuclide  vector  for  that  sample. 


Nuclide  Subsets.  In  order  to  ascertain  subsets  of  nuclides  that  are  most  effective  at  making  burnup 
predictions,  Principle  Components  Analysis  was  performed  to  see  which  nuclides  varied  the  most  in  the 
training  data  [9].  The  nuclides  that  most  contributed  to  the  first  two  principle  components  along  which 
burnup  seemed  to  vary  were  selected  and  used  to  make  predictions  on  the  test  cases.  Additionally, 
predictions  were  made  where  actinides  were  ignored. 


Effect  of  Error.  The  effect  of  measurement  error  was  investigated  by  randomly  perturbing  the  nuclide 
concentrations  in  the  test  case  nuclide  vectors  and  the  effects  on  the  accuracy  of  the  burnup  and  reactor 
type  predictions  were  noted.  For  each  element  of  a  test  case  nuclide  vector,  a  random  number  in  the 
interval  [1  -Emax,  I  +Emax]  was  chosen  to  form  an  error  vector  that  was  used  to  add  error  to  test  case  nuclide 
vectors.  The  same  error  vector  was  used  for  each  of  the  three  methods.  For  each  of  the  three  nuclide  lists 
(full,  short,  fission  products),  the  set  of  test  cases  were  analyzed  with  each  of  the  three  methods  fifty  times 
for  each  value  of  maximum  error,  i.e.  fifty  different  random  error  vectors  were  generated,  and  the  results 
averaged  in  order  to  gain  statistical  confidence  in  the  results. 


The  results  of  PCA  and  the  resulting  shortened  nuclide  list  used  in  the  analyses  are  given  in  Table  10. 
Many  of  the  identified  fission  products  reside  towards  the  edges  of  the  bimodal  fission  yield  curves. 
Additional  identified  nuclides  are  uranium,  neptunium,  and  plutonium  isotopes,  which  are  directly 
affected  by  burnup  and  neutron  flux  [13]. 
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Table  10.  Nuclide  list  derived  from  PC  A  performed  on  the  training  data. 


Nuclide 


138Ba 

103Rh 

140Ce 

101Ru 

142Ce 

102Ru 

144Ce 

104Ru 

133Cs 

150Sm 

135Cs 

88Sr 

137Cs 

90Sr 

139La 

"Tc 

95Mo 

130Te 

97Mo 

234U 

98Mo 

235pj 

100Mo 

236pj 

143Nd 

23Spj 

144Nd 

131Xe 

145Nd 

132Xe 

146Nd 

134Xe 

148Nd 

136Xe 

237Np 

89y 

iOSpd 

91Zr 

141Pr 

92Zr 

239pu 

93Zr 

240Pu 

94Zr 

241Pu 

242Pu 

96Zr 

Figure  3  shows  the  success  rate  of  the  reactor  predictions  using  all  three  nuclide  lists  and  the  nearest 
neighbor  methods.  The  results  are  averaged  over  all  the  test  cases  and  fifty  trials  each.  While  the  data  is 
extremely  noisy,  a  general  downward  trend  towards  minimum  asymptotes  at  success  rates  of 
approximately  0.2903  and  0.3548  may  be  observed.  When  data  without  added  error  is  used,  the  fission 
products  nuclide  list  exhibits  the  poorest  performance,  but  lines  corresponding  to  these  lists  seems  to 
approach  the  asymptote  slower  than  the  other  nuclide-list/method  combinations,  which  implies  that 
measurements  including  actinides  seems  to  perform  the  best  but  are  most  sensitive  to  measurement 
inaccuracies.  Using  the  two  lists  including  actinides,  in  order  to  achieve  better  than  50%  success  in 
reactor  type  predictions,  the  maximum  error  in  nuclide  measure  must  be  less  than  1%. 
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Figure  3.  The  success  rate  at  predicting  the  reactor  class  (PWR,  BWR,  PHWR)  on  the  test  cases.  Data 
shown  is  averaged  over  the  62  test  cases  and  50  trials/unique  error  vectors,  which  were  used  for  both  of 
the  reactor  type  prediction  methods  as  described  in  Section  3.2. 

Figures  4  and  5  show  the  percent  error  in  bumup  predictions  averaged  over  the  test  set  and  fifty  trials. 
For  small  amounts  of  added  error,  ridge  regression  makes  the  most  accurate  predictions  and  does  so  until 
the  maximum  measurement  error  exceeds  1  or  2  percent.  Not  including  actinides  in  the  nuclide  lists  with 
the  nearest  neighbor  approaches  improves  the  algorithms’  accuracy  with  clean  data  compared  to  the  other 
nuclide  lists.  Additionally,  the  decline  in  performance  with  increasing  measurement  error  is  slowest  using 
this  “fission  products”  nuclide  list.  Furthermore,  using  the  Li  norm  seems  to  yield  better  performance 
than  the  L2  norm  at  all  levels  of  test-case  measurement  error. 
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Figure  4.  Average  percent  error  in  burnup  predictions  on  the  test  cases.  Data  shown  is  averaged  over  the 
62  test  cases  and  50  trials/unique  error  vectors,  which  were  used  for  all  of  the  predictions  methods  as 
described  in  Section  3.2. 
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Figure  5.  Detailed  average  percent  error  in  burnup  prediction  on  the  test  cases,  zoomed  into  the 
region  between  0  and  1  percent. 


While  this  work  is  promising,  it  is  by  no  means  exhaustive.  We  outline  several  limitations  of  this  study  as 
disclaimers  and  suggestions  for  further  work. 


Analysis  Methods.  There  is  significant  optimization  that  could  be  done  with  the  methods  employed,  such 
as  extending  the  first  two  analysis  methods  described  in  Section  2  to  use  more  similar  training  data,  i.e., 
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k-nearest  neighbors  [9],  Additionally,  there  are  other  feature  selection  methods  that  could  be  employed  to 
separate  nuclides  that  are  useful  for  predictions  from  those  that  are  confusing,  such  as  correlation  analysis 
or  sequential  feature  extraction  methods,  which  could  improve  the  predictive  power  of  the  methods 
employed  in  this  study  [9].  Moreover,  there  are  numerous  approaches  to  multivariate  data  analysis  that 
offer  many  additional  benefits  to  the  simple  approaches  taken  in  this  initial  study,  such  as  noise  reduction 
and  robustness  to  outliers  and  measurement  error  [9], 

Data  Preparation.  There  has  been  no  examination  in  this  work  on  the  effect  of  using  a  smaller  training 
set  or  measurement  error  in  the  training  data,  and  how  these  issues  could  be  mitigated,  such  as  replicate 
measurements  taken  during  the  collection  of  training  data.  Also,  the  nuclide  lists  used  in  this  study  were 
formed  without  thought  about  detection  methods  and  which  radionuclides  could  be  measured  with 
reasonable  accuracy  and  precision.  Furthermore,  the  artificial  error  added  to  the  test  cases  was  done  with 
a  random  distribution.  A  non-trivial  improvement  would  be  to  add  error  sampled  from  a  Poisson 
distribution  while  taking  into  account  parent/daughter  relationships. 

Through  this  simple  study  of  nearest  neighbor  and  ridge  regression  algorithms  to  burnup  and  reactor  class 
predictions  using  simulated  nuclear  composition  vectors  of  simulated  spent  nuclear  fuel,  it  has  been 
demonstrated  that  good  results  can  be  obtained  when  the  measurements  are  accurate  to  within  1%.  These 
results  suggest  that  investigation  into  refining  these  techniques  by  applying  other  analysis  algorithms  and 
more  extensive  testing  multivariate  data  analysis’s  application  to  this  and  other  forensics  problems  is 
warranted. 

Task  4.9  Mass  spectrometry  measurement  of  146Sm  and  151Sm. 


Radio-samarium  samples  were  created  for  measurement  in  The  University  of  Texas  TRIGA 
reactor.  15  is  a  nearly  pure  beta  emitter  apart  from  one  low  energy,  low  intensity  gamma  ray. 
It  is  difficult  to  detect  this  gamma  peak  in  activated  natural  Sm  which  contains  appreciable 
Compton  interference  and  x-rays  from  the  other  activation  products.  Beta  spectroscopy  is  also 
problematic  as  self-shielding  effects  and  interferences  from  6  beta  rays  of  the  154Sm  activation 
product,  155Eu,  complicate  the  quantification  of  the  underlying  151Sm  beta  spectrum.  The 
alternative  approach  to  direct  measurement  compares  measured  values  of  a  second  activation 
product  in  the  sample  with  solutions  of  the  Bateman  equations. 


The  151Sm  content  may  be  ascertained  from  gamma-ray  spectra  of  155Eu  which  is  produced  from 
the  activation  of  154Sm  and  subsequent  decay  of  155Sm.  The  Bateman  equations  describing  this 
set  of  reactions  are 


dN 


Sm  151 


dt 


~  ^y,Sml50^Sml50  (Fa,Sml5l  +  As»i151  , 


Sml51 


dN, 


Sm  155 


dt 


=  f  m  —  (p  +;  w  =  F  N  -X  N 
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dN 


Eu  155 


dt 


~  ^Sml55^ Sml55  ^Eul55^ Eul55 


(12) 


The  N  terms  are  number  densities,  the  F  terms  are  radiative  capture  rates,  the  Fa  terms  are 

absorption  rates,  and  A  are  decay  constants.  Note  the  approximation  in  equation  12.  The  22.3 
min  half-life  of  155Sm  results  in  a  decay  constant  which  is  several  orders  of  magnitude  larger  than 
the  total  absorption  rate.  Thus  the  contribution  to  losses  due  to  absorption  is  negligible.  It  is 
assumed  that  for  short  to  moderate  irradiation  times  (days  or  less)  the  number  densities  of  150Sm 
and  154Sm  remain  very  nearly  constant.  These  equations  may  then  be  solved  using  Laplace 
transforms.  The  number  density  of 15 'Sm  is 


F  a I 

jy  (j.\  _  y,Sm\50  Sm\ 50  _  ~(^a,Sml5l+^Sml5\)t  ~l 

y  Sm  151V1/  j-,  q  Lx  J 

r a,Sm\5\  +  Asml51 


The  155Eu  number  density  is  given  by 


(13) 
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Natural  isotope  number  densities  and  decay  constants  were  obtained  from  [14].  Capture  and 
absorption  rates  were  calculated  from 

<15> 

where  i  denotes  the  particular  type  of  interaction  (e.g.  y,151Sm),  <j>  is  the  flux  in  the  energy  group 

g,  o t ,,  is  the  cross  section  for  interaction  of  type  i  in  energy  group  g,  and  G  is  the  total  number  of 

energy  groups.  63  group  cross  sections  from  the  Cinder90  cross  section  library  included  in 
MCNPX  2.6  [15]  were  used  as  the  or  .  Group  fluxes  were  calculated  by  means  of  MCNPX 

using  a  volumetric  (F4)  flux  tally  in  a  cell  containing  a  cylinder  of  Sm(III)  oxide.  The  size, 
position,  and  density  of  the  cell  were  defined  to  replicate  the  material  and  geometric  description 
of  the  sample  irradiated  in  the  TRIGA  core.  Using  equation  7  is  equivalent  to  using  1  group 
collapsed  cross  sections  in  the  activation  equations. 

The  irradiation  consists  of  the  irradiation  of  a  0.5  g  of  Sm(III)  powder  sample  was  irradiated  in  a 
rotary  specimen  rack  facility  (RSR)  of  the  University  of  Texas  (UT)  at  Austin’s  TRIGA  Mark  II 
reactor.  The  RSR  is  a  rotating  annular  channel  surrounding  the  fuel  assembly.  The  powder  was 
irradiated  for  30  minutes  at  a  reactor  power  of  950kW.  This  corresponds  to  a  total  neutron 
fluence  of  5X 1015  cm’2. 


34 


Quantification  of  the  155Eu  concentration  was  done  by  gamma  ray  spectroscopy.  After  decaying 
for  several  days,  the  irradiated  Sm  powder  was  counted  on  a  high  purity  gennanium  detector 
(HPGe)  for  24  hours  in  normal  mode.  Areas  of  the  86.5,  105.3  and  146.1  keV  photopeaks  were 
used  to  extrapolate  the  155Eu  activity  in  the  sample  at  the  end  of  irradiation.  In  order  to  determine 
the  decay  rates  of  the  86.5  and  105.3  keV  peaks,  a  low-energy  efficiency  curve  was  generated 
with  calibrated  137Cs  and  57Co  sources. 

The  151Sm  concentration  was  also  estimated  using  bumup  calculations  within  MCNPX 
2.6.  A  model  of  UT  Austin’s  TRIGA  core  was  modified  to  include  a  BURN  card.  The  BURN 
card  couples  burnup  and  activation  calculations  with  a  KCODE  criticality  calculation  to 
iteratively  adjust  spectral  changes  in  the  flux  and  track  the  buildup  of  various  nuclides.  See  [15] 
for  more  details.  To  reduce  the  computation  time  required  to  arrive  at  reasonable  uncertainties, 
variance  reduction  was  implemented  by  increasing  the  neutron  importance  of  cells  near  the  Sm 
cell.  This  includes  a  number  of  nested  virtual  cells  which  split  each  incoming  neutron  into  two. 
Color  maps  of  the  neutron  importance  are  shown,  imposed  over  the  core  geometry,  in  figure  6. 
The  virtual  cells  are  the  nested  cylinders  surrounding  the  Sm  sample  (the  smallest  blue  cell).  All 
other  cells  represent  real  geometric  elements  of  the  core. 


Figure  6:  Neutron  importance  in  the  MCNPX  model  of  the  UT  Austin’s  TRIGA  core,  (a)  side 
view  and  (b)  top  view  of  sample  in  the  rotary  specimen  rack  (RSR)  irradiation  facility. 


The  Sm  oxide  powder  was  counted  for  24  hours  after  a  29  day  decay  period.  The  activity 
concentration  of  155Eu  extrapolated  to  the  end  of  irradiation  period  was  17.7±0.2  kBq/g. 
Equations  4  and  5  were  used  to  estimate  a  151Sm/150Sm  ratio  of  (1.35±0.15)x10‘7  corresponding 
to  the  measured  concentration  of  155Eu.  The  MCNPX  burnup  calculations  estimated  the  ratio  to 
be  (1.3±0.2)xl0‘7  and  the  155Eu  concentration  to  be  20±3  kBq/g.  The  calculated  values  for  the 
155Eu  concentration  and  the  Sm  ratio  are  shown  as  the  solid  black  line  in  figure  7  and  represent 
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solutions  to  the  Bateman  equations  parameterized  by  the  neutron  fluence.  The  dotted  black  lines 
indicate  the  standard  error  propagated  from  the  uncertainty  of  the  flux.  The  yellow  dot  is  located 
at  the  intersection  of  the  measured  155Eu  concentration  and  the  parametric  curve.  The 
corresponding  value  on  the  x-axis  indicates  the  Sm  ratio  inferred  from  measurement.  The  red 
square  indicates  the  values  of  the  155Eu  concentration  and  Sm  ratio  as  determined  from  the  output 
files  of  the  MCNPX  burnup  calculations.  Standard  error  bars  on  l55Eu  concentration  and  Sm 
ratio  are  displayed  for  each  data  point  in  figure  7.  The  error  on  the  measured  1j5Eu  concentration, 
however,  is  too  small  to  display.  The  combined  average  of  both  estimates  is  (1.33±0.15)x10"7. 

The  MCNPX  model  ran  30  inactive  KCODE  cycles  and  130  active  cycles  with  10000  sources 
particles  per  cycle.  The  run  took  18  hours  on  a  Pentium  T4400  dual-core  processor. 


Figure  7:  Sm  ratio  estimates  made  from  measured  values  and  MCNPX  burnup  calculations. 


The  Sm  ratio  estimated  from  the  measured  155Eu  concentrations  is  in  reasonable  agreement  with 
the  value  calculated  with  MCNPX.  Large  uncertainties  are  associated  with  the  MCNPX  results 
despite  variance  reduction  efforts.  This  is  attributed  to  the  small  size  of  the  sample  relative  to  the 
reactor  geometry  which  results  in  a  small  probability  of  the  Sm  cell  interacting  with  a  source 
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neutron.  Consequently  the  solution  converges  slowly  relative  to  flux  tallies  in  other  cells  but  still 
orders  of  magnitude  faster  than  with  no  variance  reduction  present. 

Several  naturally  occurring  Sm  isotopes  have  unusually  large  absorption  cross  sections.  149Sm 
for  example  has  a  60232b  Maxwell  averaged  cross  section.  Therefore,  even  a  sample  of  Sm 
oxide  powder  as  small  as  1  gram  can  significantly  perturb  the  neutron  flux.  The  effects  of 
neutron  self-shielding  greatly  alter  the  spatially  and  energy  dependent  flux  from  its  values  in  an 
infinitely  dilute  absorber.  Using  the  155Eu  concentration  to  estimate  the  Sm  ratio  is  akin  to  using 
a  flux  monitor  with  the  advantage  that  the  homogeneously  distributed  154Sm  accounts  for  any 
large  spatial  and  spectral  variations  in  the  flux  that  contribute  to  the  lj0Sm  activation  rate. 
Therefore,  using  155Eu  to  renonnalize  the  flux  only  requires  an  accurate  knowledge  of  the 
spectral  shape  in  the  neutron  flux.  This  is  usually  within  the  capabilities  of  MCNPX. 

151Sm,  prepared  by  neutron  activation  of  Sm(III)oxide  (SrrnCE)  powder  was  quantified  for  future 
accelerator  mass  spectrometry  measurements.  Ratios  of  151Sm/150Sm  estimated  by  gamma 
spectrometry  measurements  of  155Eu  coupled  with  solutions  of  Bateman  equations  agreed  with 
MCNPX  burnup  calculations.  These  two  estimates  yielded  an  average  Sm  ratio  of 
(1.33±0.15)xl0‘7. 

The  irradiated  SrrnCE  powder  was  shipped  to  LLNL  for  isotopic  analyses.  It  was  dissolved  in 
Ultrex  HC1  (<  1x10" 12  Eu  and  Nd)  to  a  concentration  of  8  mg  Sm/mL  in  5%  HC1.  An  aliquot  of 
this  stock  solution  was  used  to  establish  the  detection  limits  by  ICPMS  described  in  Task  4.7. 

The  stock  solution  was  diluted  with  Sm  ICP  standard  (<1x10  6  Eu)  solution  to  produce  solutions 
of  different  151Sm/150Sm  ratios.  These  diluted  solutions,  ICP  standard  solution,  and  non- 
irradiated  Sm203  dissolved  in  HC1  were  prepared  as  SmF3  samples  as  described  in  Task  4.4  with 
slight  changes  described  below. 


Table  11.  Isotope  ratios  of  151Sm/150Sm  in  Sm  solutions. 


Sample  Name 

Sm  concentration  (mg/mL) 

151Sm/15USm 

Working  Stock 

8.1 

1.33xl0'7 

Dilution  1 

1.47 

4.88xl0'8 

Dilution  2 

1.06 

3.43xl0"9 

Dilution  3 

1.01 

4.73x10"'° 

Dilution  4 

1.01 

4.29x10" 

Dilution  5 

1.00 

3.67xl0"12 

Dilution  6 

1.00 

6.29xl0"13 

Start  with  sufficient  volume  to  have  3  mg  Sm  in  solution.  Add  3  mL  48%  ultrex  hydrofluoric 
acid  (<lxl0~12  Eu  and  Nd)  and  allow  samarium  fluoride  to  precipitate  overnight.  Centrifuge  to 
form  a  pellet  and  remove  the  supernatant.  Add  1  mL  of  dilute  hydrofluoric  acid,  re-suspend  the 
pellet  and  centrifuge  again.  Repeat  rinse  with  dilute  hydrofluoric  acid  and  centrifuge  2  more 
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times.  Dry  precipitate  overnight  in  a  convection  oven  set  at  85°C.  Place  dry  SmF3  in  a  watertight 
vial  until  ready  to  use.  The  day  before  measurement,  mix  each  SmF3  sample  with  Nb  metal 
powder  (Sm:Nb~3: 1)  and  press  into  aluminum  sample  targets.  The  mass  of  SmF3  was  measured 
and  indicated  Sm  recovery  was  75-80%.  Targets  were  also  loaded  with  commercially  available 
high  purity  SmF3  powder  (<lxl06  Eu  and  Nd)  mixed  with  Nb  and  used  to  tune  the  AMS  system 
for  the  Sm  beam.  Place  the  targets  in  the  ion  source  overnight  to  further  dry  the  samples  under 
vacuum.  The  ion  source  pressure  is  <lxl0~6  Torr. 

The  AMS  system  was  tuned  using  the  152SmF-  produced  by  the  ion  source  operating  at  a  cathode 
voltage  of  5.5  kV  and  extraction  voltage  of  -40  kV.  Approximately  80  uA  of  negative  ions  are 
produced  by  the  source  with  80  nA  being  152SmF-  after  a  low  energy  mass  spectrometer.  The 
'“Sm  isotope  was  selected  because  neither  Sm  nor  "  Sm  are  found  naturally.  Thus,  Sm  is 
well  separated  from  tails  of  other  stable  Sm  isotopes.  The  152SmF-  beam  was  accelerated 
through  a  potential  of  7.5MV  and  sent  through  a  thin  foil  to  strip  electrons  to  produce  a 
distribution  of  positive  charge  states.  At  7.5MV,  the  most  populated  charge  states  are  9+  and 
10+,  both  at  about  20%.  We  selected  152Sm9+  because  odd  numbered  charge  states  tend  to  have 
fewer  interferences,  so  our  high  energy  beam  had  an  energy  of  74.2 1  MeV.  Substantial 
scattering  occurs  during  passage  through  the  stripping  foil,  so  the  focused  beam  intensity 
between  the  accelerator  and  high  energy  mass  spectrometer  is  only  about  150  nA  despite  each 
ion  carrying  a  charge  of  9+.  After  the  high  energy  mass  spectrometer,  the  152Sm91  drops  to  only 
1-2  nA,  and  is  maintained  through  a  electrostatic  analyzer  to  the  particle  detector.  The 
electrostatic  analyzer  is  essentially  an  energy  filter.  These  stable  isotope  currents  were  similar  for 
AMS  targets  regardless  of  the  starting  source  of  Sm.  The  stable  152Sm9f  was  scaled  to  15 'Sm9 
and  checked  for  current  at  Faraday  cups  throughout  the  AMS  system.  No  current  was  measurable 
(<1  pA)  at  mass  151  after  the  high  energy  spectrometer.  The  slits  immediately  after  the  low 
energy  spectrometer  were  then  completely  wound  in  to  stop  the  beam  and  all  in  beam  Faraday 
cups  were  removed  from  the  beam  line.  No  counts  were  registered  in  the  detector.  The  slit  were 
then  opened  to  0.010  inches,  l/30th  the  aperture  used  for  tuning.  The  detector  registered  a  count 
rate  of  several  hundred  events  per  second,  indicative  of  a  stable  isobar  interference  (151Eu)  or 
large  tail  of  down  scattered  150Sm.  There  was  no  way  to  discern  what  caused  the  counts,  and 
opening  the  slits  any  further  would  damage  the  detector. 

The  high  count  rate  at  mass  15 1  is  very  likely  due  to  151Eu.  The  concentration  of  Eu  in  the 
various  spec  sheets  was  consistent  with  an  Eu/Sm  concentration  around  1  ppm.  This  level  is 
below  or  at  detection  limits  of  ICPMS  but  lOO-lOOOx  the  levels  tolerated  by  AMS.  Measuring 
15  'Sm  by  AMS  requires  a  significant  improvement  in  the  suppression  of  Eu,  from  the  current 
lxlO"6  available  to  lxlO"9  to  lxlO12.  Achieving  a  reduction  in  stable  Eu  concentration  of  4-6 
orders  of  magnitude  is  a  significant  task,  beyond  the  scope  of  this  limited  study.  The  lanthanides 
have  similar  chemistry  and  many  stable  isobars,  real  challenges  to  AMS.  Serial  ion 
chromatography  separations  is  an  approach,  but  it  would  be  very  labor  intensive.  The  current 
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detector  on  the  heavy  element  AMS  beamline  has  limited  dE/dx  capabilities  also,  so  an  improved 
detector  with  greater  resolution  to  resolve  scattered  tails  would  be  helpful.  An  improved  detector 
would  provide  better  separation,  but  chemistry  is  the  better  target  for  improvement.  Reducing 
interference  isobar  concentration  to  approximately  the  same  order  of  magnitude  as  the  target  rare 
isotope  is  needed  to  make  routine  15 *Sm  measurements. 


Task  4.10  Conduct  presentations/meetings  at  times  and  places  specified  in  the  contract 
schedule. 

Two  presentations  were  made  at  the  Methods  and  Application  in  Radioanalytical  Chemistry 
(MARC)  in  March  2012  and  one  presentation  was  made  at  the  2012  DTRA  Basic  Research 
Technical  Review. 

Task  4.11  Write  Option  Year  2  report. 

A  final  report  has  been  written  and  submitted  by  the  September  1,  2012  deadline. 


39 


References 


1.  Joint  Working  Group  of  the  American  Physical  Society  and  the  American  Association  for  the 

Advancement  of  Science,  Nuclear  Forensics:  Role,  State  of  the  Art,  and  Program  Needs  (2007). 

2.  MOODY,  K.,  HUTCHEON,  I.,  GRANT,  P.  Nuclear  Forensic  Analysis.  Taylor  &  Francis  (2005). 

3.  WEAVER,  C.,  BIEGALSKI,  S.,  BUCHHOLZ,  B.,  Assessment  of  non-traditional  isotopic  ratios  by 

mass  spectrometry  for  analysis  of  nuclear  activities.  Journal  of  Radioanalytical  and  Nuclear 
Chemistry  (2009). 

4.  WHITNEY,  S.,  BIEGALSKI,  S.,  Analytical  nuclear  fuel  cycles  from  isotopic  ratios  of  waste  products 

applicable  to  measurement  by  accelerator  mass  spectrometry.  Nuclear  Science  and  Engineering 
157,  200-209  (2007). 

5.  WILLMAN,  C„  HAKANSSON,  A.,  OSIFO,  0.,  BACHLIN,  A.,  SVARD,  S.,  Nondestructive  assay  of 

spent  nuclear  fuel  with  gamma-ray  spectroscopy.  Annals  of  Nuclear  Engergy  33,  427-438  (2006). 

6.  DIFOGGIO,  R.  Examination  of  some  misconceptions  about  near-infrared  analysis.  Applied 

Spectroscopy  49,  67-75  (1995). 

7.  OLIVIERI,  A.  Analytical  advantages  of  multivariate  data  processing.  Analytical  Chemistry  80,  5713- 

5720  (2008). 

8.  BEEBE,  K.,  PELL,  R.,  SEASHOLTZ.  Chemometrics:  A  Practical  Guide.  John  Wily  and  Sons  (1998). 

9.  HASTIE,  T.,  TIBSHIRANI,  R.,  FRIEDMAN,  J.  The  Elements  of  Statistical  Learning:  Data  Mining, 

Inference,  and  Predictions.  Springer  (2009). 

10.  ILAS,  G.,  GAULD,  I.  Analysis  of  decay  heat  measurements  for  BWR  fuel  assemblies,  Nuclear  Fuels 

and  Structural  Materials  for  the  Next  Generation  Nuclear  Reactors,  Reno,  Nevada,  June  4-8, 
2006.  Trans.  Am.  Nucl.  Soc.,  94,  385-387  (2006). 

11.  KNIEF,  R.  Nuclear  Engineering:  Theory  and  Technology  of  Commercial  Nuclear  Power.  American 

Nuclear  Society  (2008). 

12.  MURPHY,  B.  ORIGEN-ARP  Cross-Section  Libraries  for  Magnox,  Advanced  Gas-Cooled,  and 

VVER  Reactor  Designs,  Oak  Ridge  National  Laboratory  (2004). 

13.  Chart  of  Nuclides  (2012)  National  Nuclear  Data  Center,  Brookhaven  National  Laboratory. 

http://www.nndc.bnl.gov/chart/.  Accessed  12  February  2012 

14.  MCNPX  2.6.0  Manual  (2012)  MCNPX  Homepage. 

mcnpx.lanl.gov/docs/MCNPX_2. 6. 0_Manual.pdf.  Accessed  1  January  2012 


40 


Appendix 

Table  Al.  Results  of  the  analysis  of  the  various  test  cases  using  FOMi 


Different 

Enrichment 

Different 

Burnup 

Different 

Cooling 

Time 

Reactor 

Type  Guess 

Simulated 

Reactor 

Success 

Best-fit 

Burnup 

[MWd/MTU] 

Error 

i%] 

Interpolated 

Burnup 

[MWd/MTU] 

Error 

[%] 

Simulated 

Burnup 

[MWd/MTU] 

Cooling 

Time 

Guess 

Simulated 
Cooling  Time 

Success 

1 

'candu28' 

'CANDU28' 

1 

1290 

7.86% 

1552 

10.86% 

1400 

'lyears' 

'lyears' 

1 

1 

'candu37' 

'CANDU28' 

1290 

7.86% 

1522 

8.71% 

1400 

'7days' 

'lminutes' 

1 

'candu28' 

'CANDU28' 

1 

1290 

7.86% 

1547 

10.50% 

1400 

'30days' 

'7days' 

1 

'candu28' 

'CANDU28' 

1 

1290 

7.86% 

1552 

10.86% 

1400 

'30days' 

'30days' 

1 

1 

'candu37' 

'CANDU28' 

4740 

5.20% 

5318 

6.36% 

5000 

'lyears' 

'lyears' 

1 

1 

'candu37' 

'CANDU28' 

4740 

5.20% 

5301 

6.02% 

5000 

'7days' 

'lminutes' 

1 

'candu37' 

'CANDU28' 

4740 

5.20% 

5313 

6.26% 

5000 

'30days' 

'7days' 

1 

'candu37' 

'CANDU28' 

4740 

5.20% 

5318 

6.36% 

5000 

'30days' 

'30days' 

1 

1 

'candu28' 

'CANDU28' 

1 

10950 

0.45% 

11046 

0.42% 

11000 

'lyears' 

'lyears' 

1 

1 

'candu28' 

'CANDU28' 

1 

10950 

0.45% 

11047 

0.43% 

11000 

'lminutes 

'lminutes' 

1 

1 

'candu28' 

'CANDU28' 

1 

10950 

0.45% 

11047 

0.43% 

11000 

'7days' 

'7days' 

1 

1 

'candu28' 

'CANDU28' 

1 

10950 

0.45% 

11047 

0.43% 

11000 

'30days' 

'30days' 

1 

1 

'candu37' 

'CANDU28' 

6120 

0.00% 

6129 

0.15% 

6120 

'lyears' 

'2years' 

1 

1 

'candu37' 

'CANDU28' 

6120 

0.00% 

6129 

0.15% 

6120 

'lminutes 

'3minutes' 

1 

1 

'candu37' 

'CANDU28' 

6120 

0.00% 

6130 

0.16% 

6120 

'7days' 

'9days' 

1 

1 

1 

'candu37' 

'CANDU28' 

4740 

5.20% 

5318 

6.36% 

5000 

'lyears' 

'2years' 

1 

1 

1 

'candu37' 

'CANDU28' 

4740 

5.20% 

5302 

6.04% 

5000 

'7days' 

'3minutes' 

1 

1 

'candu37' 

'CANDU28' 

4740 

5.20% 

5315 

6.30% 

5000 

'30days' 

'9days' 

1 

'ce  16x16' 

'CE16' 

1 

1550 

8.82% 

1909 

12.29% 

1700 

'lyears' 

'lyears' 

1 

1 

'ce  16x16' 

'CE16' 

1 

1550 

8.82% 

1889 

11.12% 

1700 

'7days' 

'lminutes' 

1 

'ce  16x16' 

'CE16' 

1 

1550 

8.82% 

1905 

12.06% 

1700 

'30days' 

'7days' 

1 

'ce  16x16' 

'CE16' 

1 

1550 

8.82% 

1909 

12.29% 

1700 

'30days' 

'30days' 

1 

1 

'vverlOOO' 

'CE16' 

9150 

5.17% 

8307 

4.52% 

8700 

'lyears' 

'lyears' 

1 

1 

'werlOOO' 

'CE16' 

9150 

5.17% 

8317 

4.40% 

8700 

'lminutes 

'lminutes' 

1 

1 

'vverlOOO' 

'CE16' 

9150 

5.17% 

8334 

4.21% 

8700 

'lminutes 

'7days' 

1 

'vverlOOO' 

'CE16' 

9150 

5.17% 

8309 

4.49% 

8700 

'30days' 

'30days' 

1 

1 

'w  14x14' 

'CE16' 

16750 

1.47% 

17198 

1.16% 

17000 

'lyears' 

'lyears' 

1 

1 

'w  14x14' 

'CE16' 

16750 

1.47% 

17206 

1.21% 

17000 

'lminutes 

'lminutes' 

1 

1 

'w  14x14' 

'CE16' 

16750 

1.47% 

17196 

1.15% 

17000 

'30days' 

'7days' 

1 

'w  14x14' 

'CE16' 

16750 

1.47% 

17201 

1.18% 

17000 

'30days' 

'30days' 

1 

1 

1 

'vverlOOO' 

'CE16' 

9150 

5.17% 

8308 

4.51% 

8700 

'lyears' 

'2years' 

1 

1 

1 

'vverlOOO' 

'CE16' 

9150 

5.17% 

8317 

4.40% 

8700 

'lminutes 

'3minutes' 

1 

1 

1 

'vverlOOO' 

'CE16' 

9150 

5.17% 

8313 

4.45% 

8700 

'7days' 

'9days' 

1 

1 

'ce  16x16' 

'CE16' 

1 

9150 

0.00% 

9154 

0.04% 

9150 

'lyears' 

'2years' 

1 

1 

'ce  16x16' 

'CE16' 

1 

9150 

0.00% 

9164 

0.15% 

9150 

'lminutes 

'3minutes' 

1 

1 

'ce  16x16' 

'CE16' 

1 

9150 

0.00% 

9152 

0.02% 

9150 

'30days' 

'9days' 

1 

1 

'gelO' 

'CE16' 

8880 

2.07% 

8467 

2.68% 

8700 

'7days' 

'7days' 

1 

1 

1 

1 

'gelO' 

'CE16' 

8880 

2.07% 

8464 

2.71% 

8700 

'30days' 

'9days' 

1 

'gelO' 

'CE16' 

8880 

2.95% 

9125 

0.27% 

9150 

'30days' 

'7days' 

1 

1 

'gelO' 

'CE16' 

8880 

2.95% 

9126 

0.26% 

9150 

'30days' 

'9days' 

1 

'ge7' 

'GE7' 

1 

1980 

1.00% 

2062 

3.10% 

2000 

'lyears' 

'lyears' 

1 

1 

'ge7' 

'GE7' 

1 

1980 

1.00% 

2054 

2.70% 

2000 

'lminutes 

'lminutes' 

1 

1 

'ge7' 

'GE7' 

1 

1980 

1.00% 

2060 

3.00% 

2000 

'7days' 

'7days' 

1 

1 

'ge7' 

'GE7' 

1 

1980 

1.00% 

2061 

3.05% 

2000 

'30days' 

'30days' 

1 

1 

'ge8' 

'GE7' 

7500 

4.17% 

6922 

3.86% 

7200 

'lyears' 

'lyears' 

1 

1 

'ge8' 

'GE7' 

7500 

4.17% 

6913 

3.99% 

7200 

'lminutes 

'lminutes' 

1 

1 

'ge8' 

'GE7' 

7500 

4.17% 

6932 

3.72% 

7200 

'lminutes 

'7days' 

1 

'ge8' 

'GE7' 

7500 

4.17% 

6938 

3.64% 

7200 

'lminutes 

'30days' 

1 

'ge7' 

'GE7' 

1 

10950 

1.39% 

10666 

1.24% 

10800 

'lyears' 

'lyears' 

1 

1 

'ge7' 

'GE7' 

1 

10950 

1.39% 

10658 

1.31% 

10800 

'lminutes 

'lminutes' 

1 

1 

'ge7' 

'GE7' 

1 

10950 

1.39% 

10662 

1.28% 

10800 

'7days' 

'7days' 

1 

1 

'ge7' 

'GE7' 

1 

10950 

1.39% 

10669 

1.21% 

10800 

'7days' 

'30days' 

1 

'ge7' 

'GE7' 

1 

8880 

0.00% 

8910 

0.34% 

8880 

'lyears' 

'2years' 

1 

1 

'ge7' 

'GE7' 

1 

8880 

0.00% 

8900 

0.23% 

8880 

'lminutes 

'3minutes' 

1 

1 

'ge7' 

'GE7' 

1 

8880 

0.00% 

8906 

0.29% 

8880 

'7days' 

'9days' 

1 

1 

1 

'ge8' 

'GE7' 

7500 

4.17% 

6923 

3.85% 

7200 

'lyears' 

'2years' 

1 

1 

1 

'ge8' 

'GE7' 

7500 

4.17% 

6913 

3.99% 

7200 

'lminutes 

'3minutes' 

1 

1 

1 

'ge8' 

'GE7' 

7500 

4.17% 

6933 

3.71% 

7200 

'lminutes 

'9days' 

1 

'gelO' 

'GE7' 

8190 

7.77% 

8252 

7.07% 

8880 

'lminutes 

'7days' 

1 

1 

'gelO' 

'GE7' 

8190 

7.77% 

8255 

7.04% 

8880 

'lminutes 

'9days' 

1 

1 

'gelO' 

'GE7' 

6120 

15.00% 

6637 

7.82% 

7200 

'7days' 

'7days' 

1 

1 

1 

1 

'gelO' 

'GE7' 

6120 

15.00% 

6638 

7.81% 

7200 

'7days' 

'9days' 

1 
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Table  A2.  Results  of  the  analysis  of  the  various  test  cases  using  FOM2 


Different  Different 
Enrichment  Burnup 


Different 

Cooling 

Time 


Reactor 

Type  Guess 

Simulated 

Reactor 

Best-fit 

Burnup 

[MWd/MTU] 


Interpolated 

Burnup 

[MWd/MTU] 


Simulated  Cooling 
Burnup  Time 
[MWd/MTU]  Guess 


Simulated 
! Cooling  Time 
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Following  is  some  of  the  code  in  the  fonn  of  m-files  used  to  perform  Task  4.3: 


function  makeOrigenInp(File_Name,  Fuel  Type,  Enrichment,  Bumup,  ... 
S_Power,  Num_cycles,  Mod_Den,  Basis,  T_c,  T_c_units,  Out_units,... 
NRank, Location) 

%Kenneth  Dayman— University  of  Texas  at  Austin,  PNNL— May  2011 
%This  function  takes  input  parameters  and  writes  a  *.inp  file  for  input 
%into  OrigenArp  (via  Scale6)  for  batch  runs 


%Input  parameters: 

%File_Name:  Name  of  input  file  and  the  subsequent  saved  output  after 
%  running  Origen 

%Fuel_Type:  The  type  of  fuel  assembly  in  OrigenArp  shorthand, 

%  e.g.,  cel 4x1 4 

%Enrichment:  U-235  enrichment  of  the  fuel  in  percent,  e.g.,  5 
%Bumup:  Desired  final  burnup  of  the  fuel  after  N  cycles 
%S_Power:  Power  produced  per  unit  fuel  (MW/MTU)  ~  power  level  of 
%  reactor 

%Num_Cycles:  Number  of  cycles  to  use  ( 1  or  3  typically) 

%Mod_Den:  Density  of  the  moderator  (Non-variable  for  many  reactors), 
%  g/cc 

%Basis:  Amount  of  uranium  basis  in  grams 
%T_c:  Cooling  time 

%T_c_units:  Units  of  the  cooling  time  parameter  (seconds,  minutes, 

%  months,  years) 

%Out_units:  Units  of  the  plot/out,  e.g.,  grams,  Curies,  etc. 

%NRank:  Number  of  nuclides  to  be  output,  max  200 
%Location:  Directory  where  the  .inp  will  be  located  and  where  the 
%  Origen  output  will  be  saved 


%Defines  derived  variables  and  arrays  needed  in  the  .inp  file  from  the 
%input  parameters 


%Hard-Wired  Variables  for  Certain  Reactors 
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if  strcmpi(Fuel_Type,'cel4xl4') 
ModDen  =  0.7332; 
elseif  strcmpi(Fuel_Type,'ce  1 6x  1 6') 
ModDen  =  0.71; 

elseif  strcmpi(Fuel_Type,'wl4xl4') 
ModDen  =  0.7264; 
elseif  strcmpi(Fuel_Type,'sl4xl4') 
ModDen  =  0.7283; 
elseif  strcmpi(Fuel_Type,'w  15x15') 
Mod_Den  =  0.7135; 
elseif  strcmpi(Fuel_Type,'w  17x17') 
Mod_Den  =  0.723; 

elseif  strcmpi(Fuel_Type,'w  17x1 7_ofa') 
ModDen  =  0.71; 

elseif  strcmpi(Fuel_Type,'vver440(3 .6)') 
Enrichment  =  3.6; 

Mod_Den  =  0.75; 

elseif  strcmpi(Fuel_Type,'vver440(3.82)') 
Mod_Den  =  0.75; 

Enrichment  =  3.82; 

elseif  strcmpi(Fuel_Type,'vver440(4.25)') 
Mod_Den  =  0.75; 

Enrichment  =  4.25; 

elseif  strcmpi(Fuel_Type,'vver440(4.38)') 
Mod_Den  =  0.75; 

Enrichment  =  4.38; 
elseif  strcmpi(Fuel_Type,'vverlOOO’) 
Mod_Den  =  0.7145; 
elseif  strcmpi(Fuel_Type,  'agr') 

ModDen  =  1 ; 

elseif  strcmpi(Fuel_Type,  'magnox') 
ModDen  =  1 ; 

elseif  strcmpi(Fuel_Type,  'candu37') 
Mod_Den  =  0.8121; 

Enrichment  =  0.711; 
elseif  strcmpi(Fuel_Type,  'candu28') 
Mod_Den  =  0.8121; 

Enrichment  =  0.711; 
end 
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time_of_cycle  =  (Bumup/Num_cycles)/S_Power; 
power  =  S_Power*(Basis/le6); 
time_of_cycle_inc  =  time_of_cycle/10; 
tempBasis  =  le6; 
factor  =  Basis/tempBasis; 

%computes  the  number  of  atoms  of  each  U  isotope,  rounding  to  acheive  a 
%total  of  le6  atoms 

u234  =  round(((9e-05)*Enrichment+9e-9)*tempBasis); 
u235  =  round(Enrichment*tempBasis/100); 
u236  =  round(0.516854*u234); 

%changed  to  tempBasis  to  avoid  getting  a  negative  amount  of  u238 
u238  =  tempBasis  -  u236  -  u235  -  u234; 

%rescales  the  number  of  each  U  isotope  to  get  the  recover  the  desired 
%basis 

u234  =  u234*factor; 
u235  =  u23  5*  factor; 
u236  =  u23 6*  factor; 
u238  =  u23  8*  factor; 

%  Write  the  *.inp  file 

%prints  the  arp  block 

fid  =  fopen(strcat(File_Name,’.inp'),'w'); 

fprintf(fid,  -arp\r\n'); 

fprintf(fid,  ’%s',Fuel_Type); 

fprintf(fid,  W); 

fprintf(fid,  '%g', Enrichment); 

fprintf(fid,  W); 

fprintf(fid,  ’%d’,Num_cycles); 

fprintf(fid,  W); 

for  i  =  1  :Num_cycles 
if  time  of  cycle  <10 

fprintf(fid,  '%0.6f,time_of_cycle); 
fprintf(fid,'\r\n'); 
elseif  time  of  cycle  <100 

fprintf(fid,  ’%0.5f,time_of_cycle); 
fprintf(fid,'\r\n'); 
elseif  time  of  cycle  <  1000 
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fprintf(fid,  '%0.4f,time_of_cycle); 
fprintf(fid,'\r\n'); 
elseif  timeofcycle  <  10000 

fprintf(fid,  ’%0.3f,time_of_cycle); 
fprintf(fid,'\r\n'); 
end 
end 

for  i  =  1  :Num_cycles 

fprintf(fid,'%g’,S_Power); 

fprintf(fid,'\r\n'); 

end 

for  i  =  1  :Num_cycles 
fprintf(fid,  'l\r\n'); 
end 

fprintf(fid,’%g’,Mod_Den); 
fprintf(fid,'\r\n'); 
fprintf(fid,'ft33fOO  l\r\n'); 
fprintf(fid,'end\r\n'); 

%prints  the  origens  block 
fprintf(fid,’#origens\r\n'); 
fprintf(fid,'0$$  a4  33  al  1  71  e  t\r\n’); 
fprintf(fid,  ’%s',Fuel_Type); 
fprintf(fid,'\r\n'); 

fprintf(fid,’3$$  33  a3  1  27  al6  2  a33  18  e  tW); 
fprintf(fid,’35$$  0  t\r\n'); 

fprintf(fid,’56$$  1010al00al34al53al81  eW); 

fprintf(fid,’57**  0  a3  le-05  1  e\r\n'); 

fprintf(fid,’95$$  0  t\r\n'); 

fprintf(fid, ’Cycle  1 

fprintf(fid,'%s',File_Name); 

fprintf(fid,'\r\n'); 

fprintf(fid,  ’%g',Basis/le6); 

fprintf(fid,’  MTU\r\n'); 

fprintf(fid,’58**  '); 

for  i  =  1:5 

fprintf(fid,'%g’, power); 
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fprintf(fid,' '); 
end 

fprintf(fid,  '\r\n '); 
for  i  =  1:5 

fprintf(fid,’%g’,  power); 
fprintf(fid,’ '); 
end 

fprintf(fid,'\r\n'); 
fprintf(fid,’60**  '); 
j  =  0; 

for  i  =  1:10 
if  i  ==  6 

fprintf(fid,'\r\n '); 
end 

j  =  j  +  time_of_cycle_inc; 
if  j  <  10 

fprintf(fid,  '%0.6f  ,j); 
fprintf(fid,’  ’); 
elseifj  <100 

fprintf(fid,  '%0.5f  ,j); 
fprintf(fid,’ '); 
elseifj  <  1000 

fprintf(fid,  '%0.4f  ,j); 
fprintf(fid,' '); 
elseifj  <  10000 

fprintf(fid,  ’%0.3f  ,j); 
fprintf(fid,'  ’); 
end 
end 

if  T_c  ==  .1 
num_time_points  =  7; 
elseif  T_c  >  .1  &&  T_c  <=  .3 
num_time_points  =  6; 
elseif  T_c  >  .3  &&  T_c  <=  1 
num_time_points  =  7; 
elseif  T_c  >  1  &&  T_c  <=  3 
num_time_points  =  8; 
elseif  T_c  >  3  &&  T_c  <=  10 
num_time_points  =  9; 
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elseif  T_c  >10 

num_time_points  =  8; 
end 

fprintf(fid,'\r\n'); 

fprintf(fid,’66$$  al  2  a5  2  a9  2  e\r\n73$$  922340  922350  922360  922380W); 

fprintf(fid,  74**  '); 

fprintf(fid,  '%g',u234); 

fprintf(fid,  ’ '); 

fprintf(fid,  '%g',u235); 

fprintf(fid,  ’ '); 

fprintf(fid,  '%g',u236); 

fprintf(fid,  ’ '); 

fprintf(fid,  ’%g',u238); 

fprintf(fid,  '\r\n'); 

fprintf(fid,75$$  2  2  2  2\r\n’); 

fprintf(fid,’t\r\n'); 

fprintf(fid,’54$$  a8  1  all  0  e\r\n'); 
fprintf(fid,’56$$  a2 '); 
fprintf(fid,’%d’,num_time_points); 
fprintf(fid,'  a6  1  alO  10  ’); 
if  strcmpi(T_c_units, ’seconds’) 
fprintf(fid,  ’a  14  1  ’); 
elseif  strcmpi(T_c_units,  ’minutes’) 
fprintf(fid,  ’a  14  2  ’); 
elseif  strcmpi(T_c_units,  ’hours’) 
fprintf(fid,  ’a  14  3  ’); 
elseif  strcmpi(T_c_units,  ’years’) 
fprintf(fid,  ’a  14  5  ’); 
end 

fprintf(fid,  ’a  15  3  al7  2  e\r\n’); 
fprintf(fid,’57**  0  a3  le-05  e\r\n’); 
fprintf(fid,’95$$  0  t\r\n’); 
fprintf(fid, ’Cycle  1  Down  -  ’); 
fprintf(fid,’%s’,File_Name); 
fprintf(fid,’\r\n’); 
fprintf(fid,  ’%g’, Basis/1  e6); 
fprintf(fid,’  MTU\r\n’); 
fprintf(fid,’60**  ’); 
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%CHECK  THIS  BLOCK 
if  T_c  ==  .1 
num_time_points  =  7; 

fprintf(fid,  ’0.0001  0.0003  0.001  0.003  0.01  0.03  ’); 
fprintf(fid,'%g',T_c); 
elseif  T_c  >  .1  &&  T_c  <=  .3 
num_time_points  =  6; 
fprintf(fid,  ’0.001  0.003  0.01  0.03  0.1  ’); 
fprintf(fid,'%g’,T_c); 
elseif  T_c  >  .3  &&  T_c  <=  1 
num_time_points  =  7; 
fprintf(fid, ’0.001  0.003  0.01  0.03  0.1  0.3  '); 
fprintf(fid,'%g’,T_c); 
elseif  T_c  >  1  &&  T_c  <=  3 
num_time_points  =  8; 

fprintf(fid,  ’0.001  0.003  0.01  0.03  0.1  0.3  1  '); 
fprintf(fid,  ’%g’,T_c); 
elseif  T_c  >  3  &&  T_c  <=  10 
num_time_points  =  9; 

fprintf(fid,  ’0.001  0.003  0.01  0.03  0.1  0.3  1  3  ’); 
fprintf(fid,  ’%g’,T_c); 
elseif  T_c  >10 

num_time_points  =  8; 
fprintf(fid,  ’0.01  0.03  0.1  0.3  1  3  10  ’); 
fprintf(fid,  ’%g’,T_c); 
end 

%END  BLOCK 

fprintf(fid,'\r\n'); 
fprintf(fid,’61**  f0.05\r\n'); 
fprintf(fid, ’6  5  $  $\r\n') ; 

fprintf(fid, ’’’Gram-Atoms  Grams  Curies  Watts-All  Watts-Gamma\r\n'); 

fprintf(fid,’  3z  1  0  0  3z  3z  3z  6z\r\n'); 

fprintf(fid,’  3z  1  0  0  3z  3z  3z  6z\r\n'); 

fprintf(fid,'  3z  1  0  0  3z  3z  3z  6z\r\n'); 

fprintf(fid,’81$$  2  0  26  1  a7  200  e\r\n'); 

fprintf(fid,  ’82$$  '); 

for  i  =  1  :num_time_points 
fprintf(fid,  ’2  '); 
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end 


fprintf(fid,  ’e\r\n'); 

fprintf(fid,’83  *  *\r\n'); 

fprintf(fid,'  1.0000000e+07  8.0000000e+06  6.5000000e+06  5.0000000e+06 

4.0000000e+06\r\n'); 

fprintf(fid,’  3.0000000e+06  2.5000000e+06  2.0000000e+06  1.6600000e+06 

1.3300000e+06\r\n'); 

fprintf(fid,’  1.0000000e+06  8.0000000e+05  6.0000000e+05  4.0000000e+05 

3 .0000000e+05\r\n'); 

fprintf(fid,’  2.0000000e+05  1.0000000e+05  5.0000000e+04  1 .0000000e+04  eW); 

fprintf(fid, '  84  *  *\r\n') ; 

fprintf(fid,’  2.0000000e+07  6.4340000e+06  3.0000000e+06  1.8500000e+06V\n'); 

fprintf(fid,’  1.4000000e+06  9.0000000e+05  4.0000000e+05  1.0000000e+05 

1.7000000e+04W); 

fprintf(fid,'  3.0000000e+03  5.5000000e+02  1.0000000e+02  3.0000000e+01 

l.OOOOOOOe+OlW); 

fprintf(fid,’  3.0499900e+00  1.7700000e+00  1.2999900e+00  1.1299900e+00 

l.OOOOOOOe+OOW); 

fprintf(fid,’  8.0000000e-01  4.0000000e-01  3.2500000e-01  2.2500000e-01  9.9999850e-02\r\n’); 

fprintf(fid,’  5.0000000e-02  3.0000000e-02  9.9999980e-03  1.0000000e-05  e\r\n'); 

fprintf(fid,'t\r\n'); 


for  i  =  1  :num_time_points 
fprintf(fid,’56$$0  0al0  ’); 
fprintf(fid,'%g',i); 
fprintf(fid,’  e  t\r\n’); 
end 

fprintf(fid,'56$$  fD  t\r\n’); 

fprintf(fid,’end\r\n’); 

fprintf(fid,  -opusW) ; 

fprintf(fid,'LIBUNIT=33\r\n'); 

fprintf(fid,’TYPARAMS=NUCLIDES\r\n’); 

Outunits  =  upper(Outunits); 

fprintf(fid,’UNITS='); 

fprintf(fid, '%s' ,  Out_units) ; 

fprintf(fid,'\r\n'); 

fprintf(fid,’LIBTYPE=ALL\r\n'); 

fprintf(fid,'TIME='); 

Tcunits  =  upper(Tcunits); 
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fprintf(fid,'%s',T_c_units); 
fprintf(fid,'\r\n'); 
fprintf(fid,’NPOSITION='); 
for  i  =  1  :num_time_points 
fprintf(fid,'%g',i); 
fprintf(fid,' '); 
end 

fprintf(fid,’end\r\n'); 
fprintf(fid,’NRANK=('); 
fprintf(fid,'%d’,NRank); 
fprintf(fid,  ')\r\n') ; 
fprintf(fid,'end\r\n'); 
fprintf(fid,'#shell\r\n'); 
fprintf(fid,’copy  ft7 1  fOO  1 
fprintf(fid,'%s', Location); 
fprintf(fid,'W); 
fprintf(fid,'%s',File_Name); 
fprintf(fid,’.f7 1  "\r\n'); 
fprintf(fid,’del  ft7 1  fOO  l\r\n’); 
fprintf(fid,’end'); 

fclose(fid); 

end 


function  [Database,MasterNuclide,newNames]  =  makeDTRAdatabase() 
%Kenneth  Dayman  —  University  of  Texas  at  Austin  —  June  2011 

% 

%This  function  opens  each  ORIGEN  run  listed  in  the  batch  file,  and  compiles 
%the  nuclide  and  mass  information  into  a  usable  database  for  use  in  inverse 
%calculations  for  identifying  fuel  cycles  (defined  here  as  ordered  triples: 
%(reactor_type,  cooling  time,  burnup)) 

% 

%A  MasterNuclideList  is  compiled  from  all  the  nuclides  observed  in  the 
%output  of  the  ORIGEN  runs.  The  data  from  each  run  is  then  compiled  into 
%an  (N  x  1628)  matrix,  where  N  is  the  number  of  nuclides  in  the 
%MasterNuclideList  and  1628  is  the  number  of  runs.  If  a  particular  run 
%does  not  exhibit  a  nuclide  in  the  MasterNuclideList,  a  zero  will  be 
%inserted.  The  ordering  of  each  run  (column)  is  preserved  to  ensure  that 
%each  row  only  contains  data  relating  to  a  single  nuclide. 


51 


%Form  the  list  of  filenames  of  the  ORIGEN  output 

fid  =  fopen(’OrigenRuns.txt'); 

temp  =  textscan(fid,'%s  %s  %s’); 

names  =  temp  {3}; 

fclose(fid); 

newNames  =  cell(length(names),l); 

%Declares  Useful  static  variables 
NumRuns  =  length(names); 

NumNuclides  =  200; 

%Reads  in  the  nuclides  from  each  run  and  appends  each  unique  nuclide  not 
%already  in  MasterNuclideList  to  the  end  of  MasterNuclideList.  First  run’s 
%nuclide  list  is  set  equal  to  MasterNuclideList  to  start, 
longname  =  strcat(names(l),’._plotOOO.plf); 
fid  =  fopen(char(longname)); 

temp  =  textscan(fid,’%s  %*f  %*f  %*f  %*f  %*f  %*f  %*f  %*f  %f,’HeaderLines’,6); 

RunNuclides  =  temp  { 1 } ; 

fclose(fid); 

%Sets  MasterNuclideList  to  the  list  from  first  run,  then  deletes  the  last 
%one  (total)  and  decrements  the  NumbMasterNuclides 
MasterNuclide  =  RunNuclides; 

NumMasterNuclides  =  length(MasterNuclide); 

MasterNuclide(NumMasterNuclides)="; 

NumMasterNuclides  =  NumMasterNuclides- 1; 

%Goes  through  the  remaining  runs  and  searches  for  unique  nuclides, 

%appending  these  to  the  end  of  the  MasterNuclide 
for  i  =  2:NumRuns 
fprmtf(’%d’,i); 
fprintf(V); 

%opens  file  and  reads  in  the  nuclide  names 
longname  =  strcat(names(i),’._plotOOO.plf); 
fid  =  fopen(char(longname)); 
if  mod(i,4)==  1 

temp  =  textscan(fid,  ’%s  %*f  %*f  %*f  %*f  %*f  %*f  %*f  %*f  %f,'HeaderLines’,6); 
elseif  mod(i,4)==2 

temp  =  textscan(fid,  '%s  %*f  %*f  %*f  %*f  %*f  %*f  %*f  %f,'HeaderLines',6); 
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else 

temp  =  textscan(fid,  '%s  %*f  %*f  %*f  %*f  %*f  %*f  %f  ,'HeaderLines',6); 
end 

RunNuclides  =  temp  { 1 } ; 
fclose(fid); 

%compares  the  nuclides  in  the  run  to  all  nuclides  in  the  master  list 
for  j  =  1  :NumNuclides  %in  the  run 

for  k  =  1  :NumMasterNuclides  %master  list 
presentInMaster=0; 

if  strcmpi(RunNuclides(j),MasterNuclide(k)) 
presentInMaster=  1 ; 
break; 
end 
end 

if  presentInMaster==0 

MasterNuclide(NumMasterNuclides  +  1)  =  RunNuclides(j); 

NumMasterNuclides  =  NumMasterNuclides  +  1; 
end 
end 
end 

%Forming  Database 

Database  =  zeros(NumMasterNuclides,NumRuns); 

%the  1-loop  allows  for  the  loading  of  the  Database  matrix  to  group 
%together  all  reactor  types,  then  cooling  times,  then  burnup 
t  =  0; 
for  1  =  1 :4 

for  i  =  l:4:NumRuns 
t  =  t  +  1; 
fprintf('%d’,t); 
fprintf('\n'); 

newNames(t)  =  names(i); 

%opens  the  run  file  and  loads  the  nuclide  names  and  masses 
longname  =  strcat(names(i),'._plotOOO.plf); 
fid  =  fopen(char(longname)); 
if  mod(i,4)==l 

temp  =  textscan(fid,  ’%s  %*f  %*f  %*f  %*f  %*f  %*f  %*f  %*f  %f,'HeaderLines',6); 
elseif  mod(i,4)==2 
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temp  =  textscan(fid,  '%s  %*f  %*f  %*f  %*f  %*f  %*f  %*f  %f,'HeaderLines',6); 
else 

temp  =  textscan(fid,  '%s  %*f  %*f  %*f  %*f  %*f  %*f  %f  ,'HeaderLines',6); 
end 

RunNuclides  =  temp  { 1 } ; 

RunMass  =  temp  {2}; 
fclose(fid); 

forj  =  1  :NumMasterNuclides 
for  k  =  1  :NumNuclides 
if  strcmpi(RunNuclides(k),MasterNuclide(j)) 

Database(j,t)  =  RunMass(k); 
end  %implicitly  if  nuclide  in  master  isn't  in  run,  set  to  0 
end 
end 
end 
end 
end 


function  [reactor_type,cooling_time,guessBurnup,interpBurnup,smallFoM]  = 
IDfuelCycle(unknownFileName,numTimePoints, Database, MasterNuclide,names,FoMtype) 
%Kenneth  Dayman  —  University  of  Texas  —  June  2011 

% 

%This  function  takes  the  ORIGEN  ouput  of  an  unknown  fuel  cycle  (defined 
%for  this  work  as  an  ordered  triple  (reactor  type,  cooling_time,burnup) 

%and  makes  a  best  guess  against  a  database  of  the  Mass  outputs  from  a 
%number  of  known,  simulated  fuel  cycles.  The  best- fit  estimate  will 
%use  a  least-squares  metric  (squared  Euclidean  norm),  and  then 
%interpolates  an  estimate  for  the  burnup  using  the  "taxicab"  norm  as  a 
%distance  metric  from  the  samples  adjacent  the  Database  matrix  to  the 
%best-fit  fuel  cycle. 

% 

%unknownFileName  =  Name  of  unknown  fuel  cycle’s  ORIGEN  output  (w/o 
%  file  extension) 

%numTimePoints  =  number  of  time  points  in  the  unknown's  ORIGEN  output 
%Database  =  sorted  matrix  of  mass  of  each  nuclide  in  the  MasterNuclide 
%  list  for  each  of  the  fuel  cycles  in  names.  Rows  are 

%  nuclides  in  the  MasterNuclide  list  and  Columns  are  the 
%  fuel  cycles  in  names 
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%MasterNuclides  =  Listing  of  all  the  nuclides  accounted  for  in  the 
%  Database 

%names  =  all  the  fuel  cycles  accounted  for  in  the  Database 
%FoM  type  =  {’square’  =>  least  squares  Figure  of  Merit 
%  'abs'  =>  absolute  value  Figure  of  Merit 

%  ’relative'  =>  relative  absolute  value  Figure  of  Merit 

%  (possibly  ill-defined  due  to  divide  by  0's} 

%Reads  in  the  Names  and  masses  of  the  unknown  ORIGEN  output 
longname  =  strcat(unknownFileName,'._plotOOO.plf); 
fid  =  fopen(char(longname)); 
if  numTimePoints  ==  6 

temp  =  textscan(fid,  ’%s  %*f  %*f  %*f  %*f  %*f  %f,'HeaderLines’,6); 
elseif  numTimePoints  ==  7 

temp  =  textscan(fid,  ’%s  %*f  %*f  %*f  %*f  %*f  %*f  %f  ,'HeaderLines’,6); 
elseif  numTimePoints  ==  8 

temp  =  textscan(fid,  ’%s  %*f  %*f  %*f  %*f  %*f  %*f  %*f  %f,'HeaderLines',6); 
elseif  numTimePoints  ==  9 

temp  =  textscan(fid,  ’%s  %*f  %*f  %*f  %*f  %*f  %*f  %*f  %*f  %f,'HeaderLines',6); 
end 

RunNuclides  =  temp  { 1 } ; 

RunMass  =  temp  {2}; 
fclose(fid); 

NumMasterNuclides  =  length(MasterNuclide); 

NumNuclides  =  200; 

%Conforms  the  mass  numbers  to  MasterNuclides  order  and  full  listing. 

%If  a  nuclide  in  MasterNuclides  is  not  in  the  unkown,  then  a  0  is 
%listed.  If  a  nuclide  is  the  unknown  is  not  in  the  MasterNuclides 
%(very  unlikely),  that  nuclide  will  be  compared  against  a  0  external  of 
%the  database  in  the  calculating  of  the  figures  of  merit 
unknownFull  =  zeros(length(MasterNuclide),l); 

NumNotPresent  =  0; 

NotPresentMasses  =  zeros(l,l); 

forj  =  1  :NumMasterNuclides 
for  k  =  1  :NumNuclides 
if  strcmpi(RunNuclides(k),MasterNuclide(j)) 
unknownFull(j)  =  RunMass(k); 
end 
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end 

end 

%checks  for  stuff  in  RunNuclides  not  in  MasterNuclide  and  saves  any 
%nonPresentNuclides  in  NotPresentMasses  for  use  later  in  calc.  FoM’s 
for  k  =  1  :NumNuclides 
notThere  =  1 ; 

forj  =  1  :NumMasterNuclides 
if  strcmpi(RunNuclides(k),MasterNuclide(j)) 
notThere  =  0; 
break; 
end 
end 

if  notThere  ==  1 ; 

NotPresentMasses(NumNotPresent  +  1)  =  RunMass(k); 

NumNotPresent  =  NumNotPresent  +  1 ; 
end 
end 

%finds  FoM  for  each  known  fuel  cycle  in  the  Database 
FoM  =  zeros(l,length(names)); 
for  i  =  l:length(names) 

if  strcmpi(FoMtype, 'square') 

FoM(i)  =  sum((unknownFull-Databasc(:,i)).*(unknownFull-Databasc(:,i)));  %FoM 
elseif  strcmpi(FoMtype,’abs') 

FoM(i)  =  sum(abs(unknownFull-Database(:,i))); 
elseif  strcmpi(FoMtype, 'relative') 

FoM(i)  =  sum(abs(unknownFull-Database(:,i))./Database(:,i)); 
else  fprintf(’Invalid  Figure  of  Merit  type'); 
end 

if  NumNotPresent  ~=  0 
forj  =  1  :NumNotPresent 
if  strcmpi(FoMtype, 'square') 

FoM(i)  =  FoM(i)  +  (NotPresentMasses(j)A2);  %FoM 

elseif  strcmpi(FoMtype,’abs') 

FoM(i)  =  FoM(i)  +  abs(NotPresentMasses(j)); 
elseif  strcmpi(FoMtype, 'relative') 

fprintf('nuclides  not  found  in  Database  ignored.  FoM  undefined'); 
else 
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fprintf(’Invalid  Figure  of  Merit  type'); 
end 
end 
end 
end 

%Finds  Minimum  FoM  and  Declares  that  the  best-fit  fuel  cycle  (number  i) 
smallFoM  =  FoM(l); 
winner  =  1 ; 
for  i=2:length(FoM) 
if  FoM(i)<smallFoM 
smallFoM  =  FoM(i); 
winner  =  i; 
end 
end 

%Parses  name(i)  to  give  reactor  type  and  cooling  time 

tempparts  =  regexp(names(winner),'_','sphf); 

parts  =  tempparts  { 1 } ; 

reactor_type  =  parts(l); 

cooling_time  =  parts(3); 

guessBurnup  =  str2double(parts(2)); 

%Interpolates  burnup  estimate  from  the  i+1  and  i-1  fuel  cycles'  burnups 
if  (guessBurnup  ~=  600  &&  guessBurnup  ~=  17700)  ||  (guessBurnup  ~=  600  &&  guessBurnup 
~=  12330) 

temptempl  =  regexp(names(winner- 1 ),'_',' 'split'); 
temp  1  =temptemp  1  { 1 } ; 

temptemp2  =  regexp(names(winner+l ),'_', 'split'); 
temp2=temptemp2  { 1 } ; 
burnl  =  str2double(templ(2)); 
burn2  =  str2double(temp2(2)); 
dl  =  abs(FoM(winner)-FoM(winner-l)); 
d2  =  abs(FoM(winner)-FoM(winner+l)); 
interpBurnup  =  (dl/(dl+d2))*(burn2-buml)  +  burnl; 
else 

interpBurnup  =  0; 
end 
end 
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